noise exercises

Feb 13 2012
Jim Mahoney

I did both of these analytically with Mathematica to handle the algrebra and draw plots.

X = 0 or 1 = sent bit
Y = 0 or 1 = received bit
f = probability of bit flip during transmit.

P(Y|X) = probability of y given x

PyGx[f_]:=(({{1-f, f}, {f, 1-f}}))

P(X) = probability of x

Px[z_]:=(({{z}, {1-z}}))

P(Y) = probability of y

Py = PyGx[f] . Px[z]

H(Y) = entropy of y

Hy=-Total[Py * Log[2,Py]]

Plot of H(Y) vs z for a given value of f

Plot[Hy/.f->0.1,{z,0,1},PlotRange->{{0,1},{0,1}},AxesLabel->{"f", "H(Y)"}]

3D plot of H(Y), f, z .
At f=0.5, Y = 0,1 with 50/50 chance no matter what is sent. 
Since this is perfectly uncertain, it also has information 1 by our definitions ... even though it has nothing to do with the sent bits.

Plot3D[Hy,{z,0,1},{f,0,1},PlotRange->{{0,1},{0,1},{0,1}},AxesLabel->{"z","f","H"}]

P(X,Y) = probability of x and y = joint probability distribution

PyAx=Table[PyGx[f][[i,j]]*Px[z][[j]],{i,1,2},{j,1,2}]

H(Y|X) = entropy of y given x 
Note that since x is given, and the noise is the same for 0’s and 1’s, it makes sense that this is independent of z.

HyGx = Simplify[Total[-Flatten[PyAx] * Flatten[Log[2,PyGx[f]]]]]

I(X;Y) = mutual entropy = channel capacity for a given P(X).

Ixy=Hy-HyGx

Plot of capacity vs z for a given value for the noise.

Plot[Ixy/.f->0.1,{z,0,1},PlotRange->{{0,1},{0,1}},AxesLabel->{"z","I(X;Y)"}]

Plot3D[Ixy,{z,0,1},{f,0,1},
PlotRange->{{0,1},{0,1},{0,1}},
BoxRatios->{1,1,0.6},
AxesLabel->{"z","f","I"}]

Finally, a plot of the max capacity (i.e. z=0.5) as a function of the noise probability.
This is the big result : how many bits of information per bit of data can we send through
this channel as a function of the probability f of a bit flip error.

Plot[Ixy/.z->0.5,{f,0,0.5},AxesLabel->{"f","I(X;Y)"}]

X = (1,2,3) = sent symbol
Y = (a,b,c) = received bit

P(Y|X) = probability of y given x

PyGx[f_]:=(({{1-f, f/2, f/3}, {f, 1-f, 2f/3}, {0, f/2, 1-f}}))

P(X) = probability of x

Px[z_,o_]:=(({{z}, {o}, {1-z-o}}))

P(Y) = probability of y

Py = PyGx[f] . Px[z,o]

H(Y) = entropy of y

Hy=-Total[Py * Log[2,Py]]

Plot of H(Y) vs (z,o) for a given value of f

Plot3D[Hy/.f->0.2,{z,0,1},{o,0,1},
BoxRatios->{1,1,0.6},
AxesLabel->{"o", "z","H(Y)"}]

P(X,Y) = probability of x and y = joint probability distribution

PyAx=Table[PyGx[f][[i,j]]*Px[z,o][[j]],{i,1,3},{j,1,3}]

{{{(1-f) z},{((f o)/2)},{(1/3) f (1-o-z)}},{{f z},{(1-f) o},{(2/3) f (1-o-z)}},{{0},{((f o)/2)},{(1-f) (1-o-z)}}}

H(Y|X) = entropy of y given x 
  ( First remove the P(Y|X) and P(Y,X) zero terms which give 0*log(0) which Mathematica doesn’t like.)

flatNoZero[x_]:=DeleteCases[Flatten[x],0];

HyGx=Total[-flatNoZero[PyAx]*Log[2,flatNoZero[PyGx[f]]]]

I(X;Y) = mutual entropy = channel capacity for a given P(X).

Ixy=Hy-HyGx

Plot of capacity vs z for a given value for the noise.

Plot3D[Ixy/.f->0.1,{z,0,1},{o,0,1},
BoxRatios->{1,1,0.7}, PlotRange->{{0,1},{0,1},{0,1.5}},
AxesLabel->{"o", "z","I(X;Y)"}]

Finally, a plot of the max capacity (i.e. z=0.5) as a function of the noise probability.
This is the big result : how many bits of information per bit of data can we send through
this channel as a function of the probability f of a bit flip error.

In this case, Mathematica isn’t doing a great job of finding the channel max I, which happens at the (z,o) value that maximizes I. 
Seems that even when doing it numerically, it keeps heading off into ranges where things aren’t defined well.

Here I’m using their FindMaximum to look for a max at a specific value for f, and then building up a table of values.
This expression says to look at Ixy where f is 0.1, starting at (z,o)=(0.1,0.1), and then vary z,o to find the highest value.
What gets returned as a solution has the maximum in its first component.

FindMaximum[Ixy/.f->0.1, {z,0.1},{o,0.1}]

So here I do that for many values of f, saving the results in a table.

result := Table[
{ff,
FindMaximum[Ixy/.f->ff,{z,0.1},{o,0.1}][[1]]},
{ff,0.02,0.4,0.02}]

result

And then we plot the table.

ListPlot[result, AxesLabel->{"f","max I(X;Y)"}]