\documentstyle[psbox,ps,11pt,fullpage]{article}

\input epsf

\begin{document}
\header{MAS.622}{2}
\problem{1}
\part{a}

Some definitions:\\
$\omega_1$: ``Good'' {\em iddb}\\
$\omega_2$: ``Bad'' {\em iddb}\\
$x_1$: NDIR measurement\\
$x_2$: NDCR measurement\\
${\bf x}$: vector of $x_1$ and $x_2$\\
so,\\
$P(\omega_1) = .8$\\
$P(\omega_2) = .2$\\
\par
The distributions of measurements for good and bad {\em iddb}s is shown in Figure \ref{gandb}.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{p3p1pa1.eps}}
\end{center}
\caption{Distributions of {\bf x} for good (left) and bad (right) {\em iddb}s}
\label{gandb}
\end{figure}

We are looking for a discriminant between good and bad, so we will use
a Bayesian discriminant function, $g({\bf x})$.

\begin{equation}
g_i({\bf x}) = {\bf x}^t{\bf W}_i{\bf x} + {\bf w}_i^t{\bf x} + w_{i0}
\end{equation}

The constants for Mr. Beeblebrox's case are best calculated by the
following matlab code:

\begin{verbatim}
sigmag = [4 0; 0 9];
sigmab = [9 0; 0 16];
mug = [0; 0];
mub = [8; 10];

Wg = (-1/2)*inv(sigmag);
Wb = (-1/2)*inv(sigmab);
wg = inv(sigmag)*mug;
wb = inv(sigmab)*mub;
wg0 = (-1/2)*mug'*inv(sigmag)*mug - (1/2)*log(norm(sigmag))+ log(.8);
wb0 = (-1/2)*mub'*inv(sigmab)*mub - (1/2)*log(norm(sigmab))+ log(.2);
\end{verbatim}

This yields values of $W_g = \left(\begin{array}{cc}-\frac{1}{8}& 0\\0& -\frac{1}{18}\\\end{array}\right)$, $W_b = \left(\begin{array}{cc}-\frac{1}{18}& 0\\0&-\frac{1}{32}\\\end{array}\right)$, $w_g = \left(\begin{array}{c}2\\\frac{10}{9}\\\end{array}\right)$, \mbox{$w_b = \left(\begin{array}{c}0\\0\\\end{array}\right)$,} $w_{g0} = -14.8773$, and $w_{b0} = -2.9957$.
\par
This yields a discriminant surface $g_g({\bf x}) = g_b({\bf x})$ of

\begin{equation}
5x_1^2 + \frac{7x_2^2}{4} - 144x_1 - 80x_2 + 855.4752 = 0
\end{equation}

which is an equation for an elipse, as seen in Figure \ref{ell}.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{oval.eps}}
\end{center}
\caption{Decision regions for {\em iddb}s.  Blue is good, red is bad.}
\label{ell}
\end{figure}

A decision function can be evaluated, $g_g({\bf x}) - g_b({\bf x}) >
0$ to identify good {\em iddb}s.  For Mr. Beeblebrox's data, if the
left hand side of the above equation is above zero, the {\em iddb} is
good.  If it is below zero, it is bad.

\part{b}

Given these added costs for misclassification, we now want to approve
(classify as $\omega_1$) a {\em iddb} only if $\frac{p({\bf x}|\omega_1)}{p({\bf x}|\omega_2)} > \frac{(\lambda_{12}-\lambda_{22}) P(\omega_2)}{(\lambda_{21}-\lambda_{11}) P(\omega_1)}$.  In Mr. Beeblerox's case (assuming he refunds the money in cases where defective units were sold):\\
$\lambda_{11} = 0$\\
$\lambda_{12} = 25000$\\
$\lambda_{21} = 5000$\\
$\lambda_{22} = 5000$\\
$P(\omega_1) = .8$\\
$P(\omega_2) = .2$\\

so we want 

\begin{equation}
\frac{p({\bf x}|\omega_1)}{p({\bf x}|\omega_2)} > \frac{20000P(\omega_2)}{5000P(\omega_1)}
\label{dec}
\end{equation}

where

\begin{equation}
p({\bf x}|\omega_i) = \frac{1}{2\pi|\Sigma_i|^\frac{1}{2}} e^{-\frac{1}{2} ({\bf x} - {\bf \mu_i})^t\Sigma_i^{-1}({\bf x} - {\bf \mu_i})}
\end{equation}

we can take the log of both sides of Equation \ref{dec}, which becomes
$g_1({\bf x}) - g_2({\bf x}) > \log(4) = 1.3863$.  This is the same
answer as part a, but the threshold is 1.3836 instead of 0.  This
means the ellipse in which a {\em iddb} is classified as good is
smaller, in order to offset the penalty of liability cases.  The new
decision curve is shown in Figure \ref{ellb}.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{ellb.eps}}
\end{center}
\caption{Decision regions for {\em iddb}s, considering liability.  Blue is good, red is bad.}
\label{ellb}
\end{figure}

\part{c}

The percentage of {\em iddb}s sold (classified as $\omega_1$) that are bad for a given threshold will be
\begin{equation}
\int_{R_1} p(\omega_2|{\bf x}) {\rm d}{\bf x} 
\end{equation}

Numerically calculating this integral over a range of thresholds
yields the curve seen in Figure \ref{thresh}.  The irregularities are
a result of the numerical integration.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{thresh.eps}}
\end{center}
\caption{Rate of defective sold units as a function of decision threshold in Equation 2}
\label{thresh}
\end{figure}

The threshold necassary to comply with the possible legislation would
be 6.8.  Depending on the penalties proscribed in the legislation, I
would propose a full fledged risk analysis, as the added sales gained
by slightly lowering the may more than cover the penalties for
exceeding the proscribed 0.01% defective rate.  The legislation under
consideration has not yet considered penalties and such analysis is
not covered under my current contract with Mr. Beeblebrox.

\problem{2}

\part{a}

Using the following matlab code, the within class scatter was
calculated.

\begin{verbatim}
for i = 1:3
	if(i == 1)
		trainset = train1;
	else if(i==2)
		trainset = train2;
	else
		trainset = train3;
	end
	end
	[rows cols] = size(trainset);
	m = mean(trainset);
	for r = 1:rows
		tmp = trainset(r,:)-m;
		Sw = Sw+(tmp'*tmp);
	end
end
\end{verbatim}

The between class scatter was calculated with:

\begin{verbatim}
Sb = 0;
allmean = mean([mean(train1) mean(train2) mean(train3)]);
for i = 1:3
	if(i == 1)
		trainset = train1;
	else if(i==2)
		trainset = train2;
	else
		trainset = train3;
	end
	end
	[rows cols] = size(trainset);
	Sb = Sb + rows * ((mean(trainset)-allmean)'*(mean(trainset)-allmean));
	
end
\end{verbatim}

I didn't like the way you told me to calculate the W matrix, so I did
{\tt [vec val] = eig(inv(Sw)*Sb); W = [vec(:,1) vec(:,2)];} (after
checking that the first to eigenvectors in vec were the ones with the
highest eigenvalues).  I did it your way to, just to make sure we both
got the same answer.  We did.
\par
The resultant projection is in Figure \ref{scatter}.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{s3p2a.eps}}
\end{center}
\caption{Class 1 is in blue x's, class 2 in red +'s, and class 3 in green x's}
\label{thresh}
\end{figure}

\part{b}

Using the estimated parameters of the gaussian distributions a
classifier was constructed.  The decision boundaries for the
classifier are seen in Figure \ref{class}.  The error rates were
relatively low; for class 1, there were no misclassifications of the
test data, for class 2 there were 3 misclassifications (all as class
3), and for class 3 there were 5 misclassifications (3 as class 2 and
2 as class 1).  The test points are seen in Figure \ref{class} as
well.  Not all points are shown, as some are off the graph.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{class.eps}}
\end{center}
\caption{Class 1 is in green, Class 2 in blue and Class 3 in red.  The
test points appear as blue x's, black +'s and white *'s for classes 1, 2 and 3, respectively}
\label{class}
\end{figure}

\part{c}

Assuming a constant covariance, the discriminants all become linear,
as seen in Figure \ref{class2}.  The noise in the corners is due to
numerical errors.  This classifier performs marginally worse than the
classifier from part b, misclassifying 2 additional class 2 points as
class 3.  This is to be expected, as we are using a generalized
covariance.

\begin{figure}[hbtp]
\begin{center}
\mbox{\epsfysize=3in \epsfbox{class2.eps}}
\end{center}
\caption{Class 1 is in green, Class 2 in blue and Class 3 in red.  The
test points appear as blue x's, black +'s and white *'s for classes 1, 2 and 3, respectively}
\label{class2}
\end{figure}

\problem{3}

I have chosen to discuss scenario 1.
\par

If the cow photos are not all taken in the same position, such that
they are in the same place in the frame, I would apply a match filter
of the average cow over the space of translations and small rotations
to get the cow image into a standard configuration.
\par

Given this image, I would calculate the eigencows for the data set,
and project the training cow images on to this space, at which point I
would calculate a Fisher linear discriminant to distinguish the cows.
\par

Assuming all the cows have different names and are all equally likely
to be in a situation which demands their identification, and equal
prior would be assigned across all possible cows.  Of course, if some
set of the cows (say the bulls) are less likely to be brought into the
arch for identification, their prior would adjusted appropriately.
\par

Costs for the decisions would be assigned, most likely one-zero costs,
unless there were particular errors that were more critical to avoid.
\par

Testing would be done by taking a second set of test photos of the
cows and cross-validating.  This entire exercise is done under the
assumption that the cows are in fact distinguishable based on the
patterns on their backs.  It is certainly possible that many cows
would be essentially identical.  If the rancher has primarily
Holsteins, it may work out ok, but if he's got primarily Angus's, he
may be in trouble and should consider radiotagging.  Further, drastic
changes in lighting or cleanliness of the cattle could cause the
system to fail altogether.  Finally, a cow trained in the system as a
calf that matures could become unidentifiable or easily misclassified.
\par

The calf problem could be addressed by not including them in the
original sample, or adding some scaling to the initial filter,
assuming that the pattern on the cow's back will scale as the cow
does.  The lighting and cleanliness issues are most practically
addressed by improving mechanisms to assure consistent lighting, and
maybe a hose to clean the cattle.  Failing that, training data should
include as much in the way of expected variations in these parameters
as possible.
\par

If the cows are not in fact distinguishable, as it might be in the
case of Angus', other measurements might be taken, or ideally a
radiotagging scheme.  Unique brands could provide an identifiable
feature at lower cost.

\end{document}