Monday, February 28, 2022

Hyperbolic vertices

Suppose you have two distinct points on the x-axis of the coordinate plane. If I tell you a parabola passes through those two points, where on the plane could that parabola’s vertex be? Spoiler alert: The vertex can be anywhere on the perpendicular bisector of those two points. (Neat!)

Now, suppose the two distinct points are anywhere on the coordinate plane. If I tell you that a parabola with a vertical line of symmetry passes through those two points, where on the plane could that parabola’s vertex be?

Let $P(x_1, y_1)$ and $Q(x_2,y_2)$ be the two distinct points anywhere on the coordinate plane. Let's assume that $x_1 \ne x_2 \ne 0.$ Any three points in the $xy$-plane define a unique parabola with a vertical line of symmetry that passes through all three points. So, for instance assume that the third point is $R(0,t)$ for some $t \in \mathbb{R}.$ The parabola is given by $$\begin{align*}0 &= \begin{vmatrix} x^2 & x & y & 1 \\ x_1^2 & x_1 & y_1 & 1 \\ x_2^2 & x_2 & y_2 & 1 \\ 0 & 0 & t & 1 \end{vmatrix} \\ &= x^2 \left(x_1y_2 - x_2y_1 + t(x_2-x_1)\right) \\ &\quad\quad - x\left(x_1^2y_2 - x_2^2 y_1 + t(x_2^2 - x_1^2)\right) \\ &\quad\quad\quad\quad + y \left(x_1^2 x_2 - x_2^2 x_1\right) - t \left(x_1^2 x_2 - x_2^2 x_1 \right)\end{align*}$$ or equivalently $$y = \frac{x_1y_2 - x_2y_1 + t(x_2 - x_1)}{x_1x_2 (x_2 - x_1)} x^2 - \frac{x_1^2 y_2 - x_2^2 y_1 + t (x_2^2 - x_1^2)}{x_1x_2 (x_2 - x_1)} x + t.$$

The vertex of this parabola is at the point $$V(t) = V\left( \frac{x_1^2 y_2 - x_2^2 y_1 + t (x_2^2 - x_1^2) }{2 (x_1y_2 - x_2y_1 + t (x_2 - x_1))}, t - \frac{\left((x_1^2 y_2 - x_2^2 y_1 + t (x_2^2 - x_1^2) \right)^2}{4x_1x_2(x_2 - x_1) \left(x_1y_2 - x_2y_1 + t (x_2 - x_1)\right)} \right).$$ Let $x = x(t) = \frac{x_1^2 y_2 - x_2^2 y_1 + t (x_2^2 - x_1^2) }{2 (x_1y_2 - x_2y_1 + t (x_2 - x_1))},$ then we have $$y = y(t) = t - \frac{x_1y_2 - x_2y_1 + t (x_2 - x_1)}{x_1x_2(x_2-x_1)} x(t)^2.$$ By inverting $x = x(t),$ we get $$t = \frac{x_1^2 y_2 - x_2^2 y_1 - 2x(x_1y_2 - x_2y_1)}{(x_2-x_1) (2x - (x_1+x_2))},$$ so we get \begin{align*}y &= \frac{x_1^2y_2 - x_2^2y_1 - 2x(x_1y_2 - x_2y_1)}{(x_2-x_1) (2x - (x_1+x_2))} - \frac{(x_1y_2-x_2y_1) (2x - (x_1 + x_2)) - \left(x_1^2 y_2 -x_2^2 y_1 - 2x(x_1y_2 - x_2y_1)\right)}{x_1x_2(x_2 - x_1) (2x - (x_1 + x_2))} x^2 \\ &= \frac{x_1^2 y_2 - x_2^2 y_1 - 2(x_1y_2-x_2y_1) x + (y_2-y_1)x^2}{2(x_2-x_1) x - (x_2^2 - x_1^2)}\end{align*} which is a hyperbola with vertical asymptote at the midpoint between $x_1$ and $x_2$ and a second asymptote along the line $$y = \frac{1}{2} \frac{y_2-y_1}{x_2-x_1} x + \frac{(x_2y_2 - x_1y_1) -3(x_1y_2-x_2y_1)}{4(x_2-x_1)}.$$

Note that this covers the original case when $y_1 = y_2 = 0.$

Yet another misuse of the term "growing exponentially"

I have a mystery function that passes through the points $(0, 1)$, $(1, 2)$, $(2, 4)$, $(3, 8)$ and $(4, 16)$. I also know the function is continuous and smooth. In other words, you can draw it in a single stroke without any sharp corners or cusps. At this point, you might think this function has to be f(x) = 2^x. But I wouldn’t be too sure of that.

Can you find a continuous, smooth function that passes through those five points, but is decreasing at $x = 2$?

The snarky answer is of course, sure, of course I can, but if pressed I suppose I should have a good idea of how to actually do so. One way to do it is to find a polynomial $p(x)$ (which makes it automatically smooth and continuous) which satisfies $p(i) = 2^i$ for $i = 0, 1, 2, 3, 4$ and such that $p^\prime(2) \lt 0.$ For instance, to make things concrete we can have $p^\prime(2) = -1 \lt 0.$ In this case, since we have $6$ equations, we should have a $5$th degree polynomial $$p(x) = a_0 + a_1 x + a_2 x^2 + a_3 x^3 + a_4 x^4 + a_5 x^5.$$

Translating the conditions into linear equations with respect to the coefficients $a = (a_0, \dots, a_5) \in \mathbb{R}^6$ we get the system of equations \begin{align*} p(0) &= a_0 + 0 a_1 + 0 a_2 + 0 a_3 + 0 a_4 + 0 a_5 &= 1\\ p(1) &= a_0 + a_1 + a_2 + a_3 + a_4 + a_5 &= 2\\ p(2) &= a_0 + 2 a_1 + 4 a_2 + 8 a_3 + 16 a_4 + 32 a_5 &= 4\\ p(3) &= a_0 + 3 a_1 + 9 a_2 + 27 a_3 + 81 a_4 + 243 a_5 &= 8\\ p(4) &= a_0 + 4 a_1 + 16 a_2 + 64 a_3 + 256 a_4 + 1024 a_5 &= 16\\ p^\prime(2) &= 0 a_0 + a_1 + 4 a_2 + 12 a_3 + 32 a_4 + 192 a_5 &= -1\\ \end{align*} Turning this into matrix form, let $$P = \begin{pmatrix} 1 & 0 & 0 & 0 & 0 & 0 \\ 1 & 1 & 1 & 1 & 1 & 1 \\ 1 & 2 & 4 & 8 & 16 & 32 \\ 1 & 3 & 9 & 27 & 81 & 243 \\ 1 & 4 & 16 & 64 & 256 & 1024 \\ 0 & 1 & 4 & 12 & 32 & 192 \end{pmatrix} \,\, \text{ and } \,\, q = \begin{pmatrix} 1 \\ 2 \\ 4 \\ 8 \\ 16 \\ -1\end{pmatrix},$$ then the coefficients can be determined by $Pa = q.$ Since $P$ is invertible, the coefficients of the polynomial satisfying the system of equations above are given by $$a = P^{-1} q = \begin{pmatrix} 1 \\ -\frac{263}{12} \\ \frac{142}{3} \\ -\frac{1579}{48} \\ \frac{113}{12} \\ -\frac{15}{16}\end{pmatrix},$$ that is, equivalently, $$p(x) = 1 -\frac{263}{12} x + \frac{142}{3} x^2 -\frac{1579}{48} x^3 + \frac{113}{12} x^4 -\frac{15}{16} x^5.$$

Monday, February 21, 2022

Finding the luckiest coin

I have in my possession 1 million fair coins. Before you ask, these are not legal tender. Among these, I want to find the “luckiest” coin.

I first flip all 1 million coins simultaneously (I’m great at multitasking like that), discarding any coins that come up tails. I flip all the coins that come up heads a second time, and I again discard any of these coins that come up tails. I repeat this process, over and over again. If at any point I am left with one coin, I declare that to be the “luckiest” coin.

But getting to one coin is no sure thing. For example, I might find myself with two coins, flip both of them and have both come up tails. Then I would have zero coins, never having had exactly one coin.

What is the probability that I will — at some point — have exactly one “luckiest” coin?

Let $N_t$ be the number of remaining coins after $t$ flips, removing any coins that were tails. We start with $N_0 = 1,000,000.$ We could set this entire problem up to be a Markov chain with binomial transition probabilities, however, since we will want to keep track of the “luckiest” we would be better served making both $N=0$ and $N=1$ into absorbing states. That is, let $P = \{ P_{ij} \}_{i,j=1}^{N_0}$ be the transition probability matrix then $$P_{ij} = \mathbb{P} \{ N_{t+1} = j \mid N_t = i \} = \begin{cases} 1, &\text{if $i=j=0$ or $i=j=1$;}\\ \binom{i}{j} 2^{-i}, &\text{if $j \geq 2$, $0 \leq i \leq j$;}\\ 0, &\text{otherwise.}\end{cases}$$

If we define the transient transition matrix $Q$ as the portion of the transition matrix associated with states $N \in \{2, \dots, N_0\},$ then the fundamental matrix of our Markov chain is given by $N=(I-Q)^{-1}.$ Further, letting $r = (P_{21},P_{31},\dots, P_{N_01})^T = (1/2, 3/2^3, \dots, N_0/2^{N_0})^T$, then the vector $b = \{b_i\}_{i=2}^{N_0}$ given by $b = Nr$ is gives the probability of starting at node $N=i \in \{2, \dots, N_0\}$ and ending up with a single luckiest coin. This is all nice in theory, but even though this just requires solving a triangular system of linear equations, given the space and time requirements of actually calculating the answer, it is somewhat prohibitive.

Instead, we can either (a) estimate the answer using simulations or (b) rely on estimates with some sufficiently large $N_0.$ Both of these methods give the same result of roughly $72\%$ that I will find a luckiest coin.

LuckiestCoin
In [1]:
import numpy as np
from scipy.stats import binom

N0 = 1000000
def simulate(N0=N0):
    Niter = N0
    while Niter > 1:
        T = binom.rvs(Niter, 0.5, size=1)
        Niter -= T
    return Niter

num_sims = 10000
output = []
for i in range(num_sims):
    out = simulate()
    output.append(float(out))
print('{} simulations with N0={} yield an estimated probability of {}'.format(num_sims, N0, sum(output) / num_sims))


## Running matrix algebra with (N0-2) x (N0-2) matrices is dicey for N0 = 1e6, so let's
## scale it down and appeal to the fact that the answer should not change much for any N0
## as N0 -> inf
M = 5000
Q = np.eye(M-1)
for i in range(M-1):
    for j in range(i+1):
        Q[i,j] = binom.pmf(j+2, i+2, 0.5)
R = np.array([binom.pmf(1, i, 0.5) for i in range(2,M+1)])
B = np.linalg.solve(np.eye(M-1) - Q, R)
print('Using the fundamental matrix of the Markov chain gives a theoretical value of {}, for an initial number of coins of N0={}'.format(B[-1],M))
10000 simulations with N0=1000000 yield an estimated probability of 0.7181
Using the fundamental matrix of the Markov chain gives a theoretical value of 0.7213410794890949, for an initial number of coins of N0=5000
In [ ]: