Monday, November 7, 2022

Peak leaf peeping

It’s peak fall foliage season in Riddler Nation, where the trees change color in a rather particular way. Each tree independently begins changing color at a random time between the autumnal equinox and the winter solstice. Then, at a random later time for each tree — between when that tree’s leaves began changing color and the winter solstice — the leaves of that tree will all fall off at once.

At a certain time of year, the fraction of trees with changing leaves will peak. What is this maximal fraction?

Since all trees have the identical, independent peak leaf peeping distribution, in this particular case it is acceptable to lose the whole forest for a single tree. Let's define the time axis where $t = 0$ is the autumnal equinox and $t = 1$ is the winter solstice. Then we can model the randome time when this tree's leaves begins changing color as $\tau_1 \sim U(0,1).$ The second random time when the leaves all drop from the tree, $\tau_2$, is also uniformly distributed, but after $\tau_1,$ so we have $\tau_2 \sim U(\tau_1, 1).$

Here the expected probability that this tree has leaves that have changed color at time $s \in (0,1)$ (which is the same as the proportion of trees with changing leaves) is \begin{align*}P(s) = \mathbb{E} \left[ \tau_1 \leq s \lt \tau_2 \right] &= \int_0^s \int_s^1 \, \frac{du_2}{1-u_1} \, du_1 \\ &= (1-s) \int_0^s \, \frac{du_1}{1-u_1}\\ & = -(1-s) \ln (1-s).\end{align*} In order to maximize this percentage, we note that the derivative is $\frac{d}{ds} P(s) = \ln (1-s) + 1,$ so the sole critical point on $(0,1)$ is $\hat{s} = 1 - e^{-1}$ which satisfies $\frac{d}{ds} P(\hat{s}) = 0.$ At this time, the maximal proportion of trees with changing leaves is $$\hat{p} = P(1-e^{-1}) = e^{-1}.$$

Monday, March 28, 2022

Xiddlerian Algorithm

The astronomers of Planet Xiddler are back at it. They have developed a new technology for measuring the radius of a planet by analyzing its cross sections.

And so, they launch a satellite to study a newly discovered, spherical planet. The satellite sends back data about three parallel, equally spaced circular cross sections which have radii $A$, $B$ and $C$ megameters, with $0 \lt A \lt B \lt C.$ Based on these values, the scientists calculate the radius of the planet is $R$ megameters. To their astonishment, they find that $A$, $B$, $C$ and $R$ are all whole numbers!

What is the smallest possible radius of the newly discovered planet?

Let's first discuss the method to the Xiddlerian algorithmic madness. Setting up a reference frame, let's assume that we rotate the newly discovered planet so that the cross sections are parallel to the $xy$-plane and that the smallest cross section has the largest value of $z$, at say $z = k.$ Let's additionally assume that the equal distance between each cross section is $h$ for some $0 \lt h \lt R.$ Then we have the following system of equations: \begin{align*} A^2 + k^2 &= R^2 \\ B^2 + (k-h)^2 &= R^2 \\ C^2 + (k-2h)^2 &= R^2\end{align*} where $0 \lt h \lt k \lt R$ are unknown.

By taking the difference between the first two equations we get the hyperbola $$(k-h)^2 - k^2 + B^2 - A^2 = h^2 - 2kh + B^2 - A^2 = 0$$ or equivalently $$h = k \pm \sqrt{k^2 - B^2 + A^2}.$$ Naturally, since we want $h \lt k,$ we need to choose the lower branch of the hyperbola, that is $$h_1(k) = k - \sqrt{k^2 - B^2 + A^2}.$$ Similarly, by taking the difference between the first and last equations, we get the hyperbola $$(k-2h)^2 - k^2 + C^2 - A^2 = 4h^2 -4kh + C^2 - A^2 = 0$$ or equivalently $$h = \frac{k \pm \sqrt{k^2 - C^2 + A^2}}{2}.$$ Here, again, we can safely choose the lower branch since otherwise we are not guaranteed to have $(k-2h)^2 \lt (k-h)^2$ (which is equivalent to $B \lt C$), that is $$h_2(k) = \frac{k - \sqrt{k^2 - C^2 + A^2}}{2}.$$

The intersection of these two curves occurs when $$k = 2 \sqrt{k^2 - B^2 + A^2} - \sqrt{k^2 - C^2 + A^2}.$$ Squaring both sides gives $$k^2 = 4(k^2 - B^2 + A^2) - 4\sqrt{k^2 - B^2 + A^2}\sqrt{k^2 - C^2 + A^2} + (k^2 - C^2 + A^2),$$ or equivalently $$4 \sqrt{k^2 - B^2 + A^2} \sqrt{k^2 - C^2 + A^2} = 4k^2 - 4(B^2-A^2) - (C^2 - A^2).$$ Squaring both sides again gives $$16 (k^2 - B^2 + A^2) (k^2 - C^2 + A^2) = (4k^2 - 4(B^2 - A^2) - (C^2 -A^2))^2$$ or equivalently \begin{align*}0 &= 16k^4 - 16((B^2 - A^2) + (C^2 - A^2)) k^2 + (A^2 - B^2)(C^2 - A^2) \\ \quad\quad & \quad\quad\quad -16 k^4 + 8(4(B^2 - A^2) + (C^2 - A^2)) k^2 + (4(B^2 - A^2) + (C^2 - A^2))^2 \\ &= 8(2(B^2 - A^2) - (C^2 - A^2)) k^2 - (4(B^2 -A^2) + (C^2 - A^2))^2 + 16(C^2-A^2)(B^2-A^2) \\ &= 8 (2(B^2 - A^2) - (C^2 - A^2))k^2 - (4(B^2 - A^2) - (C^2 - A^2))^2.\end{align*} That is, $$k^2 = \frac{(4(B^2 - A^2) - (C^2 - A^2))^2 }{ 8(2(B^2 - A^2) - (C^2 - A^2)) } = \frac{(4B^2 -3A^2 -C^2)^2}{8(2B^2 - A^2 - C^2)}$$ and hence \begin{align*}R^2 &= A^2 + k^2 = A^2 + \frac{(4B^2 - 3A^2 - C^2)^2}{8(2B^2 - A^2 - C^2)} \\ &= \frac{8A^2 (2B^2 - A^2 - C^2) + (4B^2 - 3A^2 - C^2)^2}{8(2B^2 - A^2 - C^2)} \\ &= B^2 + \frac{(C^2 - A^2)^2}{8(2B^2 - A^2 - C^2)}\end{align*} That is, $$R = R(A,B,C) = \sqrt{ B^2 + \frac{(C^2 - A^2)^2}{8(2B^2 - A^2 - C^2)}}.$$

So we can just about start hunting and pecking for suitable triples $(A, B, C) \in \mathbb{N}^3$ such that $R \in \mathbb{N},$ but we can narrow the search by including $0 \lt A \lt B \lt C \lt R.$ Additionally, since $k^2$ must be a positive integer we must additionally have $8 \mid (C^2 - A^2)^2$ and $2B^2 - A^2 - C^2 \gt 0.$ The first conclusion is equivalent to the condition that $A \equiv C \mod 2,$ while the second conclusion ensures that we must have $C \lt \sqrt{2B^2 - A^2}$. After a while of hunting an pecking, we see that $(2, 7, 8)$ yields $R(A,B,C) = 8.$ An exhaustive search of triples $A \lt B \lt C \lt 8$ with $A \equiv C \mod 2$ and $C \lt \sqrt{2B^2 - A^2}$ reveals only four valid competitors: \begin{align*} (1,4,5) &\Rightarrow R(1, 4, 5) = 4\sqrt{7} \\ (1,6,7) &\Rightarrow R(1,6,7) = \frac{6\sqrt{165}}{11} \\ (2,5,6) &\Rightarrow R(2,5,6) = \frac{5 \sqrt{22}}{4} \\ (3,6,7) &\Rightarrow R(3,6,7) = \frac{4\sqrt{154}}{7}\end{align*} Since each of the competitors does not produce a sufficiently whole numbered $R$, the smallest radius the newly discovered planet can be is $8 = R(2,7,8)$ megameters.

Sunday, March 6, 2022

Geo and Desik's Conical Adventures

Two ants named Geo and Desik are racing along the surface of a cone. The circular base of the cone has a radius of $2$ meters and a slant height of $4$ meters. Geo and Desik both start the race on the base, a distance of 1 meter away from its center.

The race’s finish is halfway up the cone, $90$ degrees around the cone’s central axis from the start, as shown in the following diagram:

Geo and Desik both want your help in strategizing for the race. What is the length of the shortest path from the start to the finish?

Geo and Desik agree to a reference frame with which to work on devising the shortest path from start, $S,$ to finish, $F.$ They place the $x$-axis along the hyperplane that is perpendicular to the base and through the point $F,$ and the $y$-axis along the radius of the base containing $S,$ thus they conclude that $F = (-1,0,\sqrt{3})$ and $S = (0,-1,0).$

Desik decides to set up the optimization problem in two phases: (i) travel from $S$ along the circular base of the cone until you reach the conical surface and (ii) travel along the conical surface to $F.$ Desik knows that for any last point on the circular base, say $P,$ the fastest path from $S$ to $P$ is a straight line, but he is wary of assuming the same along conical surfaces.

Geo thinks geometrically and reasons that since cones have zero Gaussian curvature, he can cut along any meridian of the cone and unroll the resulting surface into a planar surface. In this case, let's say that we cut along the plane $y=0$ and then align that cut on a new $uv$-plane along the $v$-axis, such that the finish point $F$ is now located at $F=(0,2)$ and the flattened surface lies in the $u \lt 0$ halfspace. Since the slant height is the same for every point along the base of the cone, the flattened conical surface must be a circular sector of radius $R = 4.$ Furthermore, since the circumference of the base of the cone is $C = 4\pi,$ then the flattened cone must be a semi-circle of radius $R=4$ since the sector angle is $C / R = \pi.$

Geo continues explaining to Desik that in the flattened cone case, the finish point at $F = (0,2)$ on the $uv$-plane has a mirrored point $F^\prime = (0,-2),$ since in the unflatted cone case the upper and lower v-axes would be glued together. So in order to find the shortest distance from any point along the perimeter of the semi-circle, say $(-4\sin \kappa, 4\cos \kappa)$, to the $F$ we would want to take \begin{align*}d_1(\kappa) & = \min \{ \sqrt{ 16\sin^2 \kappa + (4\cos \kappa -2)^2 } , \sqrt{ 16 \sin^2 \kappa + (4 \cos \kappa + 2)^2 } \} \\ &= \min \{ 2\sqrt{5-4\cos \kappa}, 2\sqrt{5+4\cos \kappa} \}.\end{align*} The only point that Geo then needs to remind Desik of is that this angle $\kappa$ in the $uv$-plane is equivalently drawn from the apex of the cone to some point on the $xy$-plane in the original cone. In this case, since the radius of the base of the cone is $r = 2$ but the circular arc is the same length, Geo argues that in the $xy$-plane this arc should subtend an angle of $\kappa \frac{R}{r} = 2\kappa$ as in the picture below. That is, the point $(-4\sin \kappa, 4 \cos \kappa)$ in the $uv$-plane corresponds to the point $(-2\sin 2\kappa, 2\cos 2\kappa, 0)$ on the original cone. Therefore, the straight line distance from the start to this point is $$d_2(\kappa) = \sqrt{ (-2\sin 2\kappa + 1)^2 + 4\cos^2 2\kappa } = \sqrt{5 - 4\sin 2\kappa}.$$

Having successfully understood Geo's geometric handwaving, Desik sets up the optimization problem $\min \{ D(\kappa) \mid 0 \leq \kappa \lt \pi \},$ where \begin{align*}D(\kappa) &= d_1(\kappa) + d_2(\kappa)\\ &= \begin{cases} 2\sqrt{5- 4\cos \kappa} + \sqrt{5-4\sin 2\kappa}, &\text{when $\cos \kappa \geq 0$} \\ 2\sqrt{5+4\cos \kappa} + \sqrt{5-4\sin 2\kappa} , &\text{when $\cos \kappa \lt 0$.}\end{cases}\end{align*} Differentiating, we have $$\frac{dD}{d\kappa} =\begin{cases} \frac{4\sin \kappa}{\sqrt{5-4\cos \kappa}} - \frac{4 \cos 2\kappa}{\sqrt{5- 4\sin 2\kappa}}, &\text{when $\cos \kappa \gt 0$}\\ undefined, &\text{when $\cos \kappa = 0$} \\ -\frac{4 \sin \kappa}{\sqrt{5 + 4\cos \kappa}} - \frac{4 \cos 2\kappa}{\sqrt{5 -4\sin 2\kappa}}, &\text{ when $\cos \kappa \lt 0$.}\end{cases}$$ The minimum of $D$ will occur on $[0,\pi]$ either at an endpoint or a critical point, that is, where either $\frac{dD}{d\kappa} = 0$ or $\frac{dD}{d\kappa}$ is undefined. In the prior case, the only such solution within $0 \leq \kappa \lt \pi$ is $\kappa^* = \frac{\pi}{6};$ whereas, in the latter case, the only solution in $[0,\pi]$ is $\kappa^* = \frac{\pi}{2}.$ Thus evaluating $D(0) = D(\pi),$ $D(\pi/6),$ and $D(\pi/2),$ we see that Geo's and Desik's minimial distance from start to finish is $D(\pi/6) = 3\sqrt{5-2\sqrt{3}}.$

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 [ ]: