Complicating a simple Probability Problem II

Hi all, we complicated a simple problem in Probability in the previous post. Here again, we do the same thing.

The setup is this. Consider the unit interval $(0, 1)$. We start at the origin and choose $n$ random points in the interval that lies to the right. We now choose the $k$th point nearest to the origin and move there. We again choose $n$ points in the interval that lies to the right and keep moving to the right until we get inside an interval $[0, x)$ for a given $x$.

Let $\mathbb{E}_{n,k}(x)$ be the expected number of times we have to generate the numbers (or in other words, the moves we have to make). So as you may have guessed, this post about finding this number.

Okay, the procedure pretty much follows what we saw the previous time. It gets all messy and I'm not even sure how to present it. So I'm directly giving the result.

Let $P_{n,k}(z)$ denote the characteristic equation. We have,

$P_{n,k}(z)=\displaystyle\frac{z^\underline{n}}{n!}-(-1)^{n-k}\sum_{j=0}^{n-k}(-1)^j\binom{n-1-j}{k-1}\frac{z^\underline{j}}{j!}$

where $z^\underline{n}$ is the falling factorial.

Let $\eta_1,\eta_2,\cdots,\eta_n$ be the $n$ roots of the equation $P_{n,k}(z)=0$. Then we have,

$\begin{align}
\mathbb{E}_{n,k}(x)=\displaystyle\frac{H_n-\log{x}}{H_n-H_{n-k}}&+\frac{x^n}{n!}\sum_{r=1}^n\frac{x^{-\eta_r}\eta_r^\underline{n}}{\eta_rP'_{n,k}(z)}\\
&-\displaystyle\frac{x^nn^{n-1}}{n!(H_n-H_{n-k})}\sum_{r=1}^n\sum_{i=0}^{n-1}\sum_{j=0}^i\frac{(n-i)x^{-\eta_r}\eta_r^\underline{n}\eta_r^{j-1}a_{i-j}}{n^{i}P'_{n,k}(\eta_r)}
\end{align}$

where $a_m=[z^{n-m}]P_{n,k}(z)$ and $H_n$ is the $n$th Harmonic number.

Using this, we can see that

$\mathbb{E}_{n,1}(x)=1-n\log{x}$

$\mathbb{E}_{n,2}(x)=1+\displaystyle\frac{n(n-1)}{(2n-1)^2}(x^{2n-1}-1)-\frac{n(n-1)\log{x}}{2n-1}$

Considering the number of identities involving the falling factorials, there seems to enough room for simplifying the above formula but am not able to see things clearly. I think adding an animation would make it more easy to see what the 'procedure' means and am on it. I'll update it as soon as possible.


Until then
Yours Aye
Me

Complicating a simple Probability problem

This problem is all about taking a simple problem in probability and complicating it.

Consider the problem of finding the expected number of uniform variables between $0$ and $1$ needed to get a sum greater than $1$. This problem is well studied and surprisingly the answer turns out to be $e$, Euler's number.

The derivation goes something like this. Starting with $0$, we keep generating Uniform random variable $X \sim U(0,1)$ and adding them until we exceed $1$. Let $\mathbb{E}(x)$ be the expected number of random variables needed when we start from $x$.

It's easy to that $\mathbb{E}(x)=0$ for $x>1$ and $E(1)=1$, both by problem definition.

Now it's easy to write the following equation.

$\mathbb{E}(x)=1+\displaystyle\int\limits_0^1\mathbb{E}(x+t)\,dt$

Now let $y=x+t$. The only thing to take care is when $t=1$, $y=1+x$. But we already know that expected value is $0$ for a value exceeding $1$. Hence, we can only consider the integral with $1$ as the upper limit.

$\mathbb{E}(x)=1+\displaystyle\int\limits_x^1\mathbb{E}(y)\,dy$

For reasons that'll be apparent later in this post, let's define a family of functions $f_k(z)$ such that $f_k^{(k)}(z)=D^k f_k(z)=\mathbb{E}(z)$ and $f_k(1)=0$, $k > 0$. Here the bracketed exponents denotes differentiation and $D$ denotes the differentiation operator.

Using this definition, the equation above becomes,

$\begin{align}
f_1'(x)&=1+\displaystyle\int\limits_x^1\,df_1\\
&=1+f_1(1)-f_1(x)\\
&=1-f_1(x)
\end{align}$

Therefore, we end up with the differential equation $f'_1(x)+f_1(x)=1$.

Solving this with Wolfram alpha, we get $f_1(x)=1+c_1e^{-x}$.

But $\mathbb{E}(x)=f_1'(x)=-c_1 e^{-x}$. Using the boundary condition, $\mathbb{E}(1)=1$, we get $\mathbb{E}(x)=e^{1-x}$ and hence $\mathbb{E}(0)=e$.             $\blacksquare$

Now we make a slight modification. Instead of just generating one random number each time, what if we generate a bunch of $n$ uniform random numbers and choose the minimum of those to be the addend. Now we seek the expected number of random numbers added. We know that the pdf of the minimum of $n$ numbers is,

$pdf(z)=n(1-z)^{n-1}$, $0 \leq z \leq 1$

Therefore the integral equation we start with becomes,

$\mathbb{E}_{n,1}(x)=1+\displaystyle\int\limits_0^1\mathbb{E}_{n,1}(x+t)n(1-t)^{n-1}\,dt$

Note that $\mathbb{E}_{n,1}(1)=1$. Using $y=x+t$,

$\displaystyle\mathbb{E}_{n,1}(x)=1+n\int\limits_x^1\mathbb{E}_{n,1}(y)(1+x-y)^{n-1}\,dy$

Using the $f$ functions, we can transform the equation to,

$f_1'(x)=1-n f_1(x)+n(n-1)\int\limits_x^1f_1(y)(1+x-y)^{n-2}\,dy$

Now this can be repeatedly applied until the exponents inside the integral becomes $0$. We'll finally end-up with,

$\displaystyle \frac{f_n^{(n)}(x)}{n!}+\frac{f_n^{(n-1)}(x)}{(n-1)!}+\frac{f_n^{(n-2)}(x)}{(n-2)!}+\cdots+f_n(x)=\frac{1}{n!}$

Or using the differential operator,

$\displaystyle \left(\frac{D^n}{n!}+\frac{D^{n-1}}{(n-1)!}+\frac{D^{n-2}}{(n-2)!}+\cdots+1\right)f_n(x)=\frac{1}{n!}$

The characteristic equation of the differential equation is the partial exponential sum function. Let $P_{n,1}(m)$ denote this characteristic of equation and let $\eta_1,\eta_2,\cdots,\eta_n$ be the roots. That is,

$P_{n,1}(z)=\displaystyle 1+z+\frac{z^2}{2!}+\frac{z^3}{3!}+\cdots+\frac{z^n}{n!}$

Note that $P'_{n,1}(z)=P_{n-1,1}(z)$. Also, as

$P_{n,1}(z)=(z-\eta_1)(z-\eta_2)(z-\eta_3)\cdots(z-\eta_n)$

we have, $P'_{n,1}(\eta_i)=\displaystyle\prod_{i \ne j,1 \leq j \leq n}(\eta_i-\eta_j)=P_{n-1,1}(\eta_i)$

Back to the differential equation, $f_n(x)=\frac{1}{n!} + c_1 e^{\eta_1x}+ c_2 e^{\eta_2x}+ c_3 e^{\eta_3x}+\cdots+ c_n e^{\eta_nx}$

$f_n^{(k)}(x)=c_1 \eta_1^ke^{\eta_1x}+ c_2 \eta_2^ke^{\eta_2x}+ c_3\eta_3^k e^{\eta_3x}+\cdots+ c_n \eta_n^ke^{\eta_nx}$, $0\leq k\leq n$

Before solving for the constants, by this time we should realize that the constant in the RHS of the starting integral equation doesn't matter as it'll be taken care by the boundary conditions. The other constants can be solved by the conditions we enforced on the $f$ functions.

$f_n^{(k)}(1) = 0$ for $k<n$ and $f_n^{(n)}(1)=\mathbb{E}(1)=1$.

The post is getting unusually long. So I'll cut short the details. Solving for the constants using determinants, we'll have,

$\mathbb{E}_{n,1}(x)=\displaystyle\frac{1}{n!}\sum_{i=1}^n \frac{e^{-\eta_i(1-x)}\eta_i^{n-1}}{P_{n-1,1}(\eta_i)}$ and hence $\mathbb{E}_{n,1}(0)=\displaystyle\frac{1}{n!}\sum_{i=1}^n \frac{e^{-\eta_i}\eta_i^{n-1}}{P_{n-1,1}(\eta_i)}$       $\blacksquare$

Now for the final complication, what if we generate a bunch of $n$ random uniform numbers, sort them and choose the $k$th smallest number as the addend? Let $\mathbb{E}_{n,k}(x)$ be the expected value. The only thing of importance is the characteristic equation denoted as $P_{n,k}(z)$.

$P_{n,k}(z)=\displaystyle \frac{z^n}{n!}-(-1)^k \sum_{j=0}^{n-k}\binom{n-1-j}{k-1}\frac{z^j}{j!}$

That is, $P_{n,k}(D)f_n(x)=\frac{1}{n!}$

Note that still,  $P'_{n,k}(z)=P_{n-1,k}(z)$ and $P'_{n,k}(\eta_i)=\displaystyle\prod_{i \ne j,1 \leq j \leq n}(\eta_i-\eta_j)=P_{n-1,k}(\eta_i)$

And the same analysis shows,

$\mathbb{E}_{n,k}(x)=\displaystyle\frac{1}{n!}\sum_{i=1}^n \frac{e^{-\eta_i(1-x)}\eta_i^{n-1}}{P_{n-1,k}(\eta_i)}$ and hence $\mathbb{E}_{n,k}(0)=\displaystyle\frac{1}{n!}\sum_{i=1}^n \frac{e^{-\eta_i}\eta_i^{n-1}}{P_{n-1,k}(\eta_i)}$

One of the beautiful things to note here is that the nature of the final formula gives a interpretation via the contour integral. That is,

$\mathbb{E}_{n,k}(x)=\displaystyle\frac{1}{n!}\sum_{i=1}^n \frac{e^{-\eta_i(1-x)}\eta_i^{n-1}}{P_{n-1,k}(\eta_i)}=\frac{1}{2\pi i n!}\oint_{\Gamma}\frac{e^{-z(1-x)}z^{n-1}}{P_{n,k}(z)}\,dz$

where the contour $\Gamma$ is large enough to encompass all the zeroes of $P_{n,k}(z)$.

Using these I was able to find,

$\mathbb{E}_{1,1}(0)=e$, $\mathbb{E}_{2,1}(0)=e(\cos 1 + \sin 1)$ and $\mathbb{E}_{2,2}(0)=\cosh \sqrt 2$,

The other characteristic equations does not appear to have closed form solutions. But the next simpler things are $\mathbb{E}_{n,n}(0)$.

$\mathbb{E}_{n,n}(0)={}_0F_{n-1}(;1/n,2/n,\cdots,(n-1)/n;n^n/n!)$

where ${}_pF_q(z)$ is the generalized hyper-geometric function.

Solving the characteristic equations is the real trouble and I tried to find some simpler ways to do that without success. However I was able to establish,

$\mathbb{E}_{n,k}(0)\approx \displaystyle \frac{n+1}{k}\left(1+\frac{1}{2}\frac{k+1}{n+2}\right)$

This approximation is very good even when $n$ and $k$ are far apart.

I tried to cover the derivations too but the post is getting longer. Maybe I'll post them separately and break this into shorter ones.


Until then,
Yours Aye
Me

A Note on Conditional Expectations

This post is about a very small problem on conditional expectations but which got me thinking for a couple of days. Eventually I realized how I never really fully understood them.

The problems is this: Let $X, Y$ be two random variables uniformly distributed between $0$ and $1$. Find $\mathbb{E}(X\vert X+Y \leq 1)$.

I was discussing with a friend and he even tried to write the Bayes formula for exceptions. Only when I started thinking about this I realized, how may people, including me, don't see the absurdity of trying to find Expectation of events and Probability of variables.

For example, $\mathbb{E}(X+Y \leq 1)$ doesn't make any sense at all. Neither does $\mathbb{P}(X)$. Both of them completely nonsensical. For the first one to make sense, it would be more appropriate to talk about $\mathbb{E}(I_{X+Y \leq 1})$ where $I$ is an indicator random variable. Similarly, $\mathbb{P}(0.25\leq X \leq 0.5)$ is valid expression.

Let $X$ be a variable and $Z$ be an event, something like $X+Y\leq1$ or $X \geq Y$. Let $I_Z$ be the indicator random variable of the event $Z$. Now, conditioning on $Z$,

$\begin{align}\mathbb{E}(X\cdot I_z)&=\mathbb{P}(Z)\mathbb{E}(X\cdot I_z\vert Z)+\mathbb{P}(Z^c)\mathbb{E}(X\cdot I_z \vert Z^c)\\
&=\mathbb{P}(Z)\mathbb{E}(X\cdot 1\vert Z)+\mathbb{P}(Z^c)\mathbb{E}(X\cdot 0\vert Z^c)\\
&=\mathbb{P}(Z)\mathbb{E}(X\vert Z)\\
\end{align}$

On the other hand, conditioning on $X$,

$\begin{align}\mathbb{E}(X\cdot I_z)&=\mathbb{E}(\mathbb{E}(X\cdot I_z \vert X))\\
&=\mathbb{E}(X\cdot \mathbb{E}(I_z\vert X))\\
&=\mathbb{E}(X\cdot \mathbb{P}(Z\vert X))
\end{align}$

Equating the two, we finally get,

$\mathbb{E}(X\vert Z)=\displaystyle\frac{\mathbb{E}(X\cdot \mathbb{P}(Z\vert X))}{\mathbb{P}(Z)}$

Now this makes the problem tractable, as all the terms makes perfect sense provided we understand that $X$ is a variable and $Z$ is an event (which could even be something like $X=3$ in the discrete case).

An alternate way of deriving the same result will be like the following. Assume $X$ is a continuous random variable. The reasoning will be very similar if $X$ is discrete.

$\begin{align}
\mathbb{E}(X\vert Z)&=\int x f_{X\vert Z}(x)\,dx\\
&=\int x \text{ }\frac{\mathbb{P}(Z\vert X)f_{X}(x)}{\mathbb{P}(Z)}\,dx=\displaystyle\frac{1}{\mathbb{P}(Z)}\int x\text{ } \mathbb{P}(Z\vert X)f_{X}(x)\,dx\\
&=\displaystyle\frac{\mathbb{E}(X\cdot \mathbb{P}(Z\vert X))}{\mathbb{P}(Z)}
\end{align}$

However, the first way includes a lot of tiny beautiful details like the conditioning on a random variable, law of total expectation, and most importantly it gives a chance to understand that,

$\mathbb{E}(X\cdot I_z\vert I_z=1)\ne\mathbb{E}(X)$ in general. Rather when the dependence between $X$ and $Z$ is not clear or not given, it is always better to state the relation as,

$\mathbb{E}(X\cdot I_z\vert I_z=1)=\mathbb{E}(X\vert Z)$

That's what I wanted to convey in this post because I had a very difficult time on how to solve that conditional expectation. Hope you enjoyed this.


Until then
Yours Aye,
Me

Visualizing Conformal Mapping with Mathematica

Hi everyone, My first post in 2018..

Well, this post does not have any mathematics. Rather, this post is about, as the title says, visualizing conformal mapping. I recently came across an amazing Youtube Channel, 3Blue1Brown. If you haven't seen it, you should definitely check it out.

Most of the videos are pretty amazing. If you think Mathematics is beautiful with just equations, then visualizing Math is possibly the only thing better than Math itself (in my humble opinion).

Especially the video Visualizing the Riemann Zeta function is outstanding, simply the best video I've seen on Youtube, I think. Ever since I saw that, I wanted to create those animations. A simple Google search led me the fact that Grant Sanderson, creator of the Channel, uses his own Python module manim to create these animations.

I'm just a beginner in Python and using a specialized module is far from my relationship with Python. So I turned to the next best thing at my disposal, Mathematica. It took me about four evenings to get everything right to get an Animation and save it as GIF using Mathematica. But now after seeing those animations, I feel like it's really worth it.

So here's one of my first animations. The Inverse of a complex number.

Pretty neat, don't you think?

Now for the Zeta function itself. Mathematica takes about 5 mins to create the following animation.

Look at that.. So beautiful...

Moreover, looking at the two mappings we can get to see how the Zeta function has a simple pole at 1. Both graphs look almost the same at the final stage except that the Zeta function is translated by a unit distance along the real axis.

I also created one with a darker background but it doesn't look all that nice. The different colors are not by choice. In fact I'm not sure how to make all the curves along the real axis look the same color. Hope you would let me know in case you are aware of that.

Now for the code.

Clear["Global`*"];
f[z_] := Zeta[z];
Lim = 5;
DiscreteXValues =
 Join[Range[-Lim, -2], Table[j/10, {j, -15, 15}],
  Range[2, Lim]]; ContinuousYValues = {-Lim, Lim};
DiscreteYValues =
 Join[Range[-Lim, -2], Table[j/10, {j, -15, 15}],
  Range[2, Lim]]; ContinuousXValues = {-Lim, Lim};
RangeofthePlot = {{-Lim, Lim}, {-Lim, Lim}};
slides = Table[
    Show[ParametricPlot[
      Evaluate[
       Table[{XVaries (1 - t) + t Re[f[XVaries + I k]],
         k (1 - t) + t Im[f[XVaries + I k]]}, {k,
         DiscreteYValues}]], {XVaries, First[ContinuousXValues],
       Last[ContinuousXValues]}],
     ParametricPlot[
      Evaluate[
       Table[{k (1 - t) + t Re[f[k + I YVaries]],
         YVaries (1 - t) + t Im[f[k + I YVaries]]}, {k,
         DiscreteXValues}]], {YVaries, First[ContinuousYValues],
       Last[ContinuousYValues]}], PlotRange -> RangeofthePlot,
     AxesOrigin -> {0, 0}, GridLines -> Automatic], {t, 0, 1,
     1/40}]; // AbsoluteTiming
ListAnimate[slides, AnimationDirection -> ForwardBackward,
 AnimationRunning -> True, AnimationRate -> 8]
(*Export["ZetaFunction.gif",slides,ConversionRules\[Rule]{\
"ControlAppearance"\[Rule]None,"DisplayDurations"\[Rule]1/5}]*)
Export["ZetaFunction.gif",
 Flatten[{slides,
   Table[slides[[i]], {i, Length[slides] - 1, 2, -1}]}],
 "ControlAppearance" -> None, "DisplayDurations" -> 1/8]

After several iterations, I ended up with the above. It has options to control how many grid lines you want in each direction, and for each of them how long you wanna plot in the perpendicular direction. Play with DiscreteXValues, DiscreteYValues, ContinuousXValues and ContinuousYValues.

And before I finish, while we are on the 3Blue1Brown, his video on Pi hiding in prime irregularities, inspired me to find the following two 'non-obvious' results.

$\displaystyle \frac{\sqrt{2}\pi}{5}=\left\vert \frac{1}{1} + \frac{i}{2}- \frac{i}{3} - \frac{1}{4} + \frac{1}{5} + \frac{i}{6} - \frac{i}{7} - \frac{1}{8} + \cdots\right\vert$

$\displaystyle \frac{\pi}{\sqrt{10}}=\left\vert \frac{1}{1} + \frac{i}{3}- \frac{i}{7} - \frac{1}{9} + \frac{1}{11} + \frac{i}{13} - \frac{i}{17} - \frac{1}{19} + \cdots\right\vert$

Hope you enjoyed them.



Until then,
Yours Aye
Me.