%Candace Metoyer
%Stat 135 Discussion
%Chapter 8 Demos

%Example 8.10 pg 456
%95% Ellipse

clc;
clear;
format compact
format short g

load T5_8.dat
load T8_2.dat

x = T5_8;
components = T8_2;
lambda = [2770226 1429206 628129 221138 99824]';

xbar = mean(x)'
[n p] = size(x)

ReigVEC(:,1) = [0.046 0.039 -0.658 0.734 -.155]';
ReigVEC(:,2) = [-0.048 0.985 0.107 0.069 0.107]';
ReigVEC(:,3) = [0.629 -0.077 0.582 0.503 0.081]';
ReigVEC(:,4) = [-0.643 -0.151 0.250 0.397 0.586]';
ReigVEC(:,5) = [0.432 -0.007 -0.392 -0.213 0.784]';


%Plot the 95% control ellipse

alpha = 0.05;

df = 2;

chisquared = icdf('chi2', 1 - alpha, df)

%Compute axes half-lengths (following the ideas of Example 5.3)

ra = sqrt(lambda(1))*sqrt(chisquared)   %major axis  b/c lambda1 > lambda2
rb = sqrt(lambda(2))*sqrt(chisquared)   %minor axis

Axes = [ra rb];
Position = [0;0];       %Centered at (0,0)
theta12 = pi/2;         %The components are orthogonal so the angle between them is 90 degrees or pi/2.

figure
hold on
ellipse(Axes, Position, theta12, 'Major', 'off', 'Minor', 'off')
text(0,0,'+')           %Plots a plus sign at the center (0,0)


%Plot the points from the first two principal components in 2-D space with
%y2 on the vertical axis and y1 on the horizontal axis.

%Use eqn 8-22 in the text.

i = 1;
for j = 1:n
    yhat(j,i) = ReigVEC(:,i)'*(x(j,:)'-xbar);
end

i = 2;
for j = 1:n
    yhat(j,i) = ReigVEC(:,i)'*(x(j,:)'-xbar);
end

%plot the points
%make sure that your ellipse figure is still open
%to just plot the points with no ellipse, simply close the ellipse figure

plot(yhat(:,1),yhat(:,2),'.', 'MarkerSize', 15)
text(0,-3700,'{}_{^{\fontsize{10}y_1}}^{\^}')
text(-6000,0,'{}_{^{\fontsize{10}y_2}}^{\^}')
title('The 95% control ellipse based on the first two principal components')


%Construct Q-Q plots for the first two principal components

figure
qqplot(yhat(:,1))
title('Q-Q plot for first principal component')

figure
qqplot(yhat(:,2))
title('Q-Q plot for second principal component')