\underline{Problem Set 1}

MAS 622J/1.126J: Pattern Recognition and Analysis

Due Thursday, 23 September 2010

[Note: All instructions to plot data or write a program should be carried out using Matlab. In order to maintain a reasonable level of consistency and simplicity we ask that you do not use other software tools.]

If you collaborated with other members of the class, please write their names at the end of the assignment.

Problem 1: Why? [5 points]

Limit your answer to problem 1 to a page.

a. Describe an application of pattern recognition related to your research. What are the features? What is the decision to be made? Speculate on how one might solve the problem.

b. In the same way, describe an application of pattern recognition you would be interested in pursuing for fun in your life outside of work.

Solution: Refer to examples discussed in lecture.

Problem 2: Probability Warm-Up [20 points]

Let X and Y be discrete random variables. Let μ_{X} denote the expected value of X and σ_{X}^{2} denote the variance of X. a and b are constant values. Use excruciating detail to answer the following:

a. Show E[aX + bY] = aE[X] + bE[Y].

b. Show σ_{X}^{2} = E[X^{2}] − μ_{X}^{2}.

c. Show that independent implies uncorrelated.

d. Show that uncorrelated does not imply independent.

e. Let Z = aX + bY. Show that if X and Y are uncorrelated, then σ_{Z}^{2} = a^{2}σ_{X}^{2} + b^{2}σ_{Y}^{2} . f. Let X_{i} (i = 1, ..., n) be random variables independently drawn from the

same probability distribution with mean μ_{X} and variance σ_{X}^{2}. For the

sample mean \overline{X} = \frac{1}{n} ∑_{i=1}^{n} X_{i}, show the following: (i) E[\overline{X}] = μ_{X}. (ii)

Var[\overline{X}] (variance of the sample mean) = σ_{X}^{2}/n. Note that this is different

from the sample variance s_{n}^{2} = \frac{1}{n} ∑_{i=1}^{n} (X_{i} − \overline{X})^{2}.

g. Let X_{1} and X_{2} be independent and identically distributed (i.i.d) continuous random variables. Can Pr[X_{1} ≤ X_{2}] be calculated? If so, find its

value. If not, explain. Hint 1: Remember that for a continuous variable Pr[X_{1} = k] = 0, for any value of k. Hint 2: Remember the definition of i.i.d. variables.

h. Let X_{1} and X_{2} be independent and identically distributed discrete random variables. Can Pr[X_{1} ≤ X_{2}] be calculated? If so, find its value. If not, explain.

Solution:

a. The following is for continuous random variables. A similar argument holds for continuous random variables.

E[aX + bY] = ∑_{xy}(ax + by) p(x, y)

= a ∑_{xy} x p(x, y) + b ∑_{xy} y p(x, y)

= a ∑_{x} x p(x) dx + b ∑_{y} y p(y) dy
= aE[X] + bE[Y]

b. Making use of the definition of variance and the previous part, we have:

σ_{X}^{2} = E[(X − μ_{X})^{2}]
= E[X^{2} − 2μ_{X}X + μ_{X}^{2}]
= E[X^{2}] − E[2μ_{X}X] + E[μ_{X}^{2}]
= E[X^{2}] − 2μ_{X}E[X] + μ_{X}^{2}
= E[X^{2}] − 2μ_{X}μ_{X} + μ_{X}^{2}
= E[X^{2}] − 2μ_{X}^{2} + μ_{X}^{2}
= E[X^{2}] − μ_{X}^{2}

c. In order to check if two discrete random variables X and Y are uncorrelated, we have to proof σ_{xy} = 0 (the same holds for continuous random variables.)

From the previous question:

σ_{X,Y}^{2} = E[XY] − μ_{X}μ_{Y}

If two variables are independent:

E[XY] = ∑_{xy} xy p(x, y)
= ∑_{xy} xy p(x) p(y)
= ∑_{x} x p(x) ∑_{y} y p(y)
= E[X] E[Y]

Finally,

σ_{X,Y}^{2} = E[X] E[Y] − μ_{X}μ_{Y} = 0

d. To proof this, we need to find one case where 1) p(xy) ≠ p(x)p(y) and 2) σ_{xy} = 0 are satisfied. One possible solution is as follows:

Suppose we have the discrete random variables X and Y, and we observed all possibilities:

\begin{matrix}x & y \\ 1 & 1 \\ 1 & 1 \\ -1 & 1 \\ -1 & 1 \\ 0 & 0 \\ 0 & 0 \end{matrix}

If we look at the case where x = 1 and y = 1, 1) is satisfied:

p(x = 1 y = 1) = \frac{2}{6} = \frac{1}{3}
p(x = 1)p(y = 1) = \frac{2}{6} \frac{4}{6} = \frac{2}{9}

Now it is easy to verify that 2) is also satified:

σ_{X,Y}^{2} = E[XY] − μ_{X}μ_{Y} = 0 − 0\frac{4}{6} = 0

e. Given that Z = aX + bY and that X and Y are uncorrelated, we have

σ_{Z}^{2} = E[(Z − μ_{Z})^{2}]
= E[Z^{2}] − μ_{Z}^{2}
= E[(aX + bY)^{2}] − (aμ_{X} + bμ_{Y})^{2}
= E[a^{2}X^{2} + 2abXY + b^{2}Y^{2}] − (a^{2}μ_{X}^{2} + 2abμ_{X}μ_{Y} + b^{2}μ_{Y}^{2} )

= a^{2}E[X^{2}] + 2abE[XY] + b^{2}E[Y^{2}] − a^{2}μ_{X}^{2} − 2abμ_{X}μ_{Y} − b^{2}μ_{Y}^{2}
= a^{2}(E[X^{2}] − μ_{X}^{2}) + 2ab(E[XY] − μ_{X}μ_{Y}) + b^{2}(E[Y^{2}] − μ_{Y}^{2} )
= a^{2}σ_{X}^{2} + 2abσ_{XY}^{2} + b^{2}σ_{Y}^{2}
= a^{2}σ_{X}^{2} + b^{2}σ_{Y}^{2} ,

where only the last equality depends on X and Y being uncorrelated.

f. Using the result of (a) and the fact E[X_{i}] = μ_{X},

E[X̄] = E[\frac{1}{n} ∑_{i=1}^{n} X_{i}] = \frac{1}{n} ∑_{i=1}^{n} E[X_{i}] = \frac{1}{n} n μ_{X} = μ_{X}

Also, using the result of (d) and the fact Var[X_{i}] = σ_{X}^{2}

Var[X̄] = Var[\frac{1}{n} ∑_{i=1}^{n} X_{i}] = \frac{1}{n^{2}} ∑_{i=1}^{n} Var[X_{i}] = \frac{1}{n^{2}} n σ_{X}^{2} = σ_{X}^{2}/n

g. Given that X_{1} and X_{2} are continuous random variables, we know that Pr[X_{1} = x] = 0 and Pr[X_{2} = x] = 0 for any value of x. Thus,

Pr[X_{1} ≤ X_{2}] = Pr[X_{1} < X_{2}].

Given that X_{1} and X_{2} are i.i.d., we know that replacing X_{1} with X_{2} and X_{2} with X_{1} will have no effect on the world. In particular, we know that

Pr[X_{1} < X_{2}] = Pr[X_{2} < X_{1}].

However, since probabilities must sum to one, we have

Pr[X_{1} < X_{2}] + Pr[X_{2} < X_{1}] = 1.

Thus,

Pr[X_{1} ≤ X_{2}] = \frac{1}{2}.

h. For discrete random variables, unlike the continuous case above, we need to know the distributions of X_{1} and X_{2} in order to find Pr[X_{1} = x] and Pr[X_{2} = x]. Thus, the argument we used above fails. In general, it is not possible to find Pr[X_{1} ≤ X_{2}] without knowledge of the distributions of both X_{1} and X_{2}.

Problem 3: Teatime with Gauss and Bayes [20 points]

Let p(x, y) = \frac{1}{2παβ} e^{−(\frac{(y−μ)^{2}}{2α^{2}} + \frac{(x−y)^{2}}{2β^{2}} )}. a. Find p(x), p(y), p(x|y), and p(y|x). In addition, give a brief description

of each of these distributions.

b. Let μ = 0, α = 15, and β = 3. Plot p(y) and p(y|x = 9) for a reasonable range of y. What is the difference between these two distributions?

Solution:

a. To find p(y), simply factor p(x, y) and then integrate over x:

p(y) = ∫_{−∞}^{∞}p(x, y) dx

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{(y−μ)^{2}}{2α^{2}} + \frac{(x−y)^{2}}{2β^{2}} )} dx
= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{− \frac{(y−μ)^{2}}{2α^{2}}} e^{− \frac{(x−y)^{2}}{2β^{2}}} dx
= \frac{1}{\sqrt{2πα^{2}}} e^{− \frac{(y−μ)^{2}}{2α^{2}}} ∫_{−∞}^{∞} \frac{1}{\sqrt{2πβ^{2}}} e^{− \frac{(x−y)^{2}}{2β^{2}}} dx
= \frac{1}{\sqrt{2πα^{2}}} e^{− \frac{(y−μ)^{2}}{2α^{2}}}
= N (μ, α^{2})

The integral goes to 1 because it is of the form of a probability distribution integrated over the entire domain. To find p(x|y), divide p(x, y) by p(y):

p(x|y) = \frac{p(x, y)}{p(y)}
= \frac{1}{\sqrt{2πβ^{2}}} e^{− \frac{(x−y)^{2}}{2β^{2}}}
= N (y, β^{2})

Finding p(x) and p(y|x) follows essentially the same procedure, but the algebra is more involved and requires completing the square in the exponent.

p(x) = ∫_{−∞}^{∞}p(x, y) dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{(y−μ)^{2}}{2α^{2}} + \frac{(x−y)^{2}}{2β^{2}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{β^{2}(y−μ)^{2}+α^{2}(x−y)^{2}}{2α^{2}β^{2}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{β^{2}y^{2}−2β^{2}μy+β^{2}μ^{2}+α^{2}x^{2}−2α^{2}xy+α^{2}y^{2}}{2α^{2}β^{2}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{(α^{2}+β^{2})y^{2}−2(α^{2}x+β^{2}μ)y+(β^{2}μ^{2}+α^{2}x^{2})}{2α^{2}β^{2}})} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(\frac{y^{2}−2 \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}} y+ \frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(^{\frac{y^{2}−2 \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}} y+(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}−(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}+ \frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(^{\frac{(y− \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}−(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}+ \frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} dy

= ∫_{−∞}^{∞} \frac{1}{2παβ} e^{−(^{\frac{(y− \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} e^{−(^{\frac{\frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}} −(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} dy

= \frac{1}{2παβ} \sqrt{2π \frac{α^{2}β^{2}}{α^{2} + β^{2}}} e^{−(^{\frac{\frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}} −(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} ∫_{−∞}^{∞} \frac{1}{\sqrt{2π \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} e^{−(^{\frac{(y− \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )} dy

= \frac{1}{\sqrt{2π(α^{2} + β^{2})}}e^{−(^{\frac{\frac{β^{2}μ^{2}+α^{2}x^{2}}{α^{2}+β^{2}} −(\frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )}

= \frac{1}{\sqrt{2π(α^{2} + β^{2})}}e^{−(\frac{(α^{2}+β^{2})(β^{2}μ^{2}+α^{2}x^{2})−(α^{2}x+β^{2}μ)^{2}}{2α^{2}β^{2}(α^{2}+β^{2})} )}

= \frac{1}{\sqrt{2π(α^{2} + β^{2})}}e^{−(\frac{α^{2}β^{2}μ^{2}+α^{4}x^{2}+β^{4}μ^{2}+α^{2}β^{2}x^{2}−α^{4}x^{2}−2α^{2}β^{2}μx−β^{4}μ^{2}}{2α^{2}β^{2}(α^{2}+β^{2})} )}

= \frac{1}{\sqrt{2π(α^{2} + β^{2})}}e^{−(\frac{α^{2}β^{2}x^{2}−2α^{2}β^{2}μx+α^{2}β^{2}μ^{2}}{2α^{2}β^{2}(α^{2}+β^{2})} )}

= \frac{1}{\sqrt{2π(α^{2} + β^{2})}}e^{−(\frac{(x−μ)^{2}}{2(α^{2}+β^{2})})}
= N (μ, α^{2} + β^{2})

To find p(y|x) we simply divide p(x, y) by p(x). In finding p(x), we already know the form of p(y|x) (see the longest line in the derivation of p(x) above):

[width=bb=20 20 575 575]gb.png

Figure 1: The marginal p.d.f. of y and the p.d.f. of y given x for a specific value of x. Notice how knowing x makes your knowledge of y more certain.

p(y|x) = \frac{p(x, y)}{p(x)}

= \frac{1}{\sqrt{2π \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} e^{−(^{\frac{(y− \frac{α^{2}x+β^{2}μ}{α^{2}+β^{2}})^{2}}{2 \frac{α^{2}β^{2}}{α^{2}+β^{2}}}} )}

= N (\frac{α^{2}x + β^{2}μ}{α^{2} + β^{2}}, \frac{α^{2}β^{2}}{α^{2} + β^{2}} )

Note that all the above distibutions are Gaussian.

b. The following Matlab code produced Figure 1:

m = 0.0
a = 15.0
b = 3
x = 9
y = -100:1:100
mean = ((a^2)*x + (b^2)*m)/(a^2 + b^2)
var = ((a*b)^2)/(a^2 + b^2)
p_y_given_x = (1.0/sqrt(2*pi*var))*exp(-((y-mean).^2)/(2*var))
var2 = a^2 + b^2
p_y = (1.0/sqrt(2*pi*var2))*exp(-((y-m).^2)/(2*var2))
hold off
plot(y,p_y_given_x,'b')
hold on
plot(y,p_y,'r')
legend('p(y|x)', 'p(y)')
sy = size(y)
axis([y(1),y(sy(2)),0,0.2])
xlabel('y')
text(-70,0.14, '\mu=0')
text(-70,0.12, '\alpha=15')
text(-70,0.1, '\beta=3')

Problem 4: Covariance Matrix [15 points]

Let Σ = [ \begin{matrix}5 & 4 \\ 4 & 5 \end{matrix} ].

a. Find the eigenvalues and eigenvectors of Σ by hand (include all calculations.) Verify your computations with MATLAB function eig.

b. Verify that Σ is a valid covariance matrix.

c. We provide 500 data points sampled from the distribution N ([0 0], Σ). Download the dataset from the course website and project the data onto the eigenvectors of the covariance matrix. What is the effect of this projection? Include the plots before and after the projection.

Solution:

a. We can find the eigenvectors and eigenvalues of Σ by starting with the definition of an eigenvector. Namely, a vector e is an eigenvector of Σ if it satisfies

Σe = λe

for some constant scalar λ, which is called the eigenvalue corresponding to e. This can be rewritten as

(Σ − λI)e = 0.

This is equivalent to

det(Σ − λI) = 0.

Thus, we require that

(5 − λ)^{2} − 4^{2} = 0

By inspection, this is true when λ = 9 and λ = 1.

To find the eigenvectors, we plug the eigenvalues back into the equation above to get

(Σ − 9I)e = [ \begin{matrix}5 − 9 & 4 \\ 4 & 5 − 9 \end{matrix} ] [ \begin{matrix}a \\ b \end{matrix} ] = [ \begin{matrix}−4 & 4 \\ 4 & −4 \end{matrix} ] [ \begin{matrix}a \\ b \end{matrix} ] = [ \begin{matrix}0 \\ 0 \end{matrix} ] ,

which gives a = b. Normalized, this results in the eigenvector

e_{1} = [ \frac{1}{\sqrt{2}} ; \frac{1}{\sqrt{2}} ] .

Similarly, λ = 1 gives

(Σ − 1I)e = [ \begin{matrix}5 − 1 & 4 \\ 4 & 5 − 1 \end{matrix} ] [ \begin{matrix}a \\ b \end{matrix} ] = [ \begin{matrix}4 & 4 \\ 4 & 4 \end{matrix} ] [ \begin{matrix}a \\ b \end{matrix} ] = [ \begin{matrix}0 \\ 0 \end{matrix} ] ,

which gives a = −b. Normalized, this results in the eigenvector

e_{1} = [− \frac{1}{\sqrt{2}} ; \frac{1}{\sqrt{2}} ] .

b. The matrix Σ is a valid covariance matrix if it is symmetric and positive semi-definite. Clearly, it is symmetric, since Σ^{T} = Σ. One way to prove it is positive semi-definite is to show that all its eigenvalues are non-negative. This is indeed the case, as shown in the previous question.

c. Figure ?? shows how projecting the data onto the eigenvectors of the data’s covariance matrix reduces the correlation between X and Y.

Problem 5: Probabilistic Modeling [20 points]

Let x ∈ {0, 1} denote a person’s affective state (x = 0 for “positive-feeling state”, and x = 1 for “negative-feeling state”). The person feels positive with probability θ_{1}. Suppose that an affect-tagging system (or a robot) recognizes her feeling state and reports the observed state (variable y) to you. But this system is unreliable and obtains the correct result with probability θ_{2}.

a. Represent the joint probability distribution P (x, y|θ) for all x, y (a 2x2 matrix) as a function of the parameters θ = (θ_{1}, θ_{2}).

b. The Maximum Likelihood estimation criterion for the parameter θ is defined as:

θ̂_{ML} = arg max_{θ} L(t_{1}, ..., t_{n}; θ) = arg max_{θ} ∏_{i=1}^{n} p(t_{i}|θ)

where we have assumed that each data point t_{i} is drawn independently from the same distribution so that the likelihood of the data is L(t_{1}, ..., t_{n}; θ) =

∏_{i=1}^{n} p(t_{i}|θ). Likelihood is viewed as a function of the parameters, which de-

pends on the data. Since the above expression can be technically challenging, we maximize the log-likelihood log L(t_{1}, ..., t_{n}; θ) instead of likelihood. Note that any monotonically increasing function (i.e., log function) of the likelihood has the same maxima. Thus,

θ̂_{ML} = arg max_{θ} log L(t_{1}, ..., t_{n}; θ) = arg max_{θ} ∑_{i=1}^{n} log p(t_{i}|θ)

Suppose we get the following joint observations t = (x, y).

\begin{matrix}x & y \\ 1 & 0 \\ 1 & 1 \\ 0 & 0 \\ 1 & 1 \\ 0 & 0 \\ 0 & 1 \\ 0 & 0 \\ 1 & 1 \end{matrix}

What are the maximum-likelihood (ML) values of θ_{1} and θ_{2}? (Hint. Since P (x, y|θ) = P (y|x, θ_{2})P (x|θ_{1}), the estimation of the two parameters can be done separately in the log-likelihood criterion.)

Solution:

a. The probability mass function (pmf) of x ∈ {0, 1} is

P (x) = \{ \begin{matrix}θ_{1},  x = 0 \\ 1 − θ_{1},  x = 1 \end{matrix} \}

The conditional pmf of y ∈ {0, 1} given that x = 0 is

P (y|x = 0) = \{ \begin{matrix}θ_{2}, y = 0 \\ 1 − θ_{2}, y = 1 \end{matrix} \}

The conditional pmf of y given that x = 1 is

P (y|x = 1) = \{ \begin{matrix}1 − θ_{2}, y = 0 \\ θ_{2}, y = 1 \end{matrix} \}

Use P (x, y) = P (y|x)P (x) to tabulate the joint pmf of (x, y).

P (x, y) = ( \begin{matrix}P (0, 0) & P (0, 1) \\ P (1, 0) & P (1, 1) \end{matrix} ) = ( \begin{matrix}θ_{2}θ_{1} & (1 − θ_{2})θ_{1} \\ (1 − θ_{2})(1 − θ_{1}) & θ_{2}(1 − θ_{1}) \end{matrix} )

b. We select (θ_{1}, θ_{2}) to maximize the log-likelihood of the samples {(x_{i}, y_{i}), i = 1, ..., n} which may be expressed as

J(θ_{1}, θ_{2}) = ∑_{i} log P (x_{i}, y_{i})
= ∑_{i} (log P (y_{i}|x_{i}) + log P (x_{i}))

= (∑_{i} log P (y_{i}|x_{i})) + (∑_{i} log P (x_{i}))
= J_{2}(θ_{2}) + J_{1}(θ_{1})

Hence, we choose θ_{1} to maximize

J_{1}(θ_{1}) = ∑_{i} log P (x_{i})
= N(x = 1) log(1 − θ_{1}) + (n − N(x = 1)) log θ_{1}

where N(x = 1) = ∑_{i} x_{i}. Differentiating w.r.t. θ_{1} gives

\frac{∂J_{1}}{∂θ_{1}} = \frac{−N(x = 1)}{1 − θ_{1}} + \frac{n − N(x = 1)}{θ_{1}}

We set this derivative to zero and solve for θ_{1} to obtain

θ̂_{1} = 1 − \frac{N(x = 1)}{n}

Similarly, we choose θ_{2} to maximize

J_{2}(θ_{2}) = ∑_{i} log P (y_{i}|x_{i})
= N(x = y) log θ_{2} + (n − N(x = y)) log(1 − θ_{2})

where N(x = y) = ∑_{i} (x_{i}y_{i} + (1 − x_{i})(1 − y_{i})). Differentiating J_{2} w.r.t. θ_{2}, setting to zero and solving for θ_{2} gives

θ̂_{2} = \frac{N(x = y)}{n}

For the example data, θ̂_{1} = \frac{4}{8} = \frac{1}{2}, θ̂_{2} = \frac{6}{8} = \frac{3}{4}. Thus,

P̂(x, y) = ( \begin{matrix}θ̂_{2}θ̂_{1} & (1 − θ̂_{2})θ̂_{1} \\ (1 − θ̂_{2})(1 − θ̂_{1}) & θ̂_{2}(1 − θ̂_{1}) \end{matrix} )

The maximum likelihood of the data under this model is

∏_{i} P̂(x_{i}, y_{i}) = (\frac{1}{2})^{8} (\frac{1}{4})^{4} (\frac{3}{4})^{4} ≈ 4.828 × 10^{−6}

Problem 6: Ring Problem [20 points]

To get credit for this problem, you must not only write your own correct solution, but also write a computer simulation of the process of playing this game:

Suppose I hide the ring of power in one of three identical boxes while you weren’t looking. The other two boxes remain empty. After hiding he ring of power, I ask you to guess which box it’s in. I know which box it’s in and, after you’ve made your guess, I deliberately open the lid of an empty box, which is one of the two boxes you did not choose. Thus, the ring of power is either in the box you chose or the remaining closed box you did not choose. Once you have made your initial choice and I’ve revealed to you an empty box, I then give you the opportunity to change your mind – you can either stick with your original choice, or choose the unopened box. You get to keep the contents of whichever box you finally decide upon.

• What choice should you make in order to maximize your chances of receiving the ring of power? Justify your answer using Bayes’ rule.

• Write a simulation. There are two choices in this game for the contestant in this game: (1) choice of box, (2) choice of whether or not to switch. In your simulation, first let the host choose a random box to place the ring of power. Show a trace of your program’s output for a single game play, as well as a cumulative probability of winning for 1000 rounds of the two policies (1) to choose a random box and then switch and (2) to choose a random box and not switch.

Solution:

• Always switch your answer to the box you didn’t choose the first time. This reason is as follows. You have a 1/3 chance of initially picking the correct box. That is, there is a 2/3 chance the correct answer is one of the other two boxes. Learning which of the two other boxes is empty does not change these probabilities; your initial choice still has a 1/3 chance of being correct. That is, there is a 2/3 chance the remaining box is the correct answer. Therefore you should change your choice.

More formally,

event right first choice = the event that your first choice is right

event wrong first choice = the event that your first choice is wrong

event right when change = the event that you get the ring when changing your initial choice

event an empty box opened = the event that an empty box is opened after your first choice

First, P(event right first choice) = 1/3

Second, P(event right when change | event an empty box opened)

= P(event right when change, event right first choice | event an empty box opened) + P(event right when change, event wrong first choice | event an empty box opened)

= P(event right when change | event right first choice, event an empty box opened) P(event right first choice) + P(event right when change | event wrong first choice, event an empty box opened) P(event wrong first choice)

= 0*1/3 + 1*2/3 = 2/3

Thus, P(event right when change | event an empty box opened) > P(event right first choice)

Using Bayes:

C = Door with Car. O = Door Opened. S = Door Selected.

P (C|O, S) = \frac{P (O|C, S)P (C|S)}{P (O|S)}  (1)

Where:

P (O|S) = Σ_{C=1}^{3}P (O, C|S) = Σ_{C=1}^{3}P (O|C, S)P (C|S)  (2)

P (C|O, S) = \frac{1 ∗ \frac{1}{3}}{\frac{1}{2} ∗ \frac{1}{3} + 1 ∗ \frac{1}{3} + 0 ∗ \frac{1}{3}} = \frac{2}{3}  (3)

Another way to understand the problem is to extend it to 100 boxes, only one of which has the ring of power. After you make your initial choice, I then open 98 of the 99 remaining boxes and show you that they are empty. Clearly, with very high probability the ring of power resides in the one remaining box you did not initially choose.

• Here is a sample simulation output for the Ring problem:

a c t u a l : 1
guess1 : 2
r e v e a l : 3
swap : 0
guess2 : 2
a c t u a l : 3
guess1 : 3
r e v e a l : 1
swap : 0
guess2 : 3
a c t u a l : 2
guess1 : 3
r e v e a l : 1
swap : 0
guess2 : 3
swap  : 0
win  : 292
l o s e  : 708

win /( win+l o s e ) : 0.292
a c t u a l : 3
guess1 : 1
r e v e a l : 2
swap : 1
guess2 : 3
a c t u a l : 1
guess1 : 1
r e v e a l : 2
swap : 1
guess2 : 3
a c t u a l : 3
guess1 : 2
r e v e a l : 1
swap : 1
guess2 : 3
swap  : 1
win  : 686
l o s e  : 314
win /( win+l o s e ) : 0.686

Here is a Matlab program that simulates the Ring simulation output above:

f o r swap = 0 : 1
win = 0 ;
l o s e = 0 ;
f o r i = 1:1000
a c t u a l = f l o o r ( rand ()∗3)+1;
guess1 = f l o o r ( rand ()∗3)+1;
i f guess1 = a c t u a l
r e v e a l = f l o o r ( rand ()∗2)+1;
i f r e v e a l = a c t u a l
r e v e a l = r e v e a l + 1 ;

end e l s e

i f guess1 = 1 && a c t u a l = 2
r e v e a l = 3 ;
e l s e i f guess1 = 1 && a c t u a l = 3
r e v e a l = 2 ;
e l s e i f guess1 = 2 && a c t u a l = 1
r e v e a l = 3 ;
e l s e i f guess1 = 2 && a c t u a l = 3

r e v e a l = 1 ;
e l s e i f guess1 = 3 && a c t u a l = 1
r e v e a l = 2 ;
e l s e i f guess1 = 3 && a c t u a l = 2
r e v e a l = 1 ;

end end

i f swap = 1
i f guess1 = 1 && r e v e a l = 2
guess2 = 3 ;
e l s e i f guess1 = 1 && r e v e a l = 3
guess2 = 2 ;
e l s e i f guess1 = 2 && r e v e a l = 1
guess2 = 3 ;
e l s e i f guess1 = 2 && r e v e a l = 3
guess2 = 1 ;
e l s e i f guess1 = 3 && r e v e a l = 1
guess2 = 2 ;
e l s e i f guess1 = 3 && r e v e a l = 2
guess2 = 1 ;

end e l s e

guess2 = guess1 ;

end

i f guess2 = a c t u a l
win = win + 1 ;

e l s e

l o s e = l o s e + 1 ;

end

%% only p r i n t t r a c e f o r f i r s t 3 games
i f i <= 3

a c t u a l guess1 r e v e a l swap guess2

end end

%% p r i n t r e s u l t s f o r each game play p o l i c y

swap

win / ( win + l o s e )

end