clear?all
close?all
clc

%?*?Read?data
train_data?=?Creatmat（10000）;

%?Suppose?train_data?has?been?read?or?created?which?has?m?rows?and?n
%?columns.??Each?row?is?a?sample?in?which?has?n?variables/sensors.
[m?n]?=?size（train_data）;

%?*?Build?PCA?model?with?train_data
S?=?corrcoef（zscore（train_data））;

%?Standardize?train_data
pca_mean?=?mean（train_data）;
pca_std?=?std（train_data）;
[pc?score?latent]?=?princomp（zscore（train_data））;
sqrt_latent?=?sqrt（latent）;

f_t?=?zeros（1?n）;
f_q?=?zeros（1?n）;
for?i?=?1?:?n?
????%?P?is?the?loading?matrix?P?has?the?size?of?（n?x?i）
????P?=?pc（:?1?:?i）;
????%?tilde_P?has?the?size?of?（n?x?（n?-?i））
????tilde_P?=?pc（:?i?+?1?:?end）;

????%?*?Calculate?the?control?limit?of?Q?and?T^2
????%?delta_q?and?delta_t?are?control?limits?of?Q?and?T^2?indices
????beta?=?0.99;????
????delta_t?=?i?*?（m?-?1）?/?（?m?-?i）?*?finv（beta?i?m-i）;
????theta?=?zeros（3?1）;
????for?x?=?1?:?3
????????for?y?=?i?+?1?:?n
????????????theta（x）?=?theta（x）?+?latent（y）^（x）;
????????end
????end
????h?=?1?-?2?*?theta（1）?*?theta（3）?/?3?/?（theta（2）^2）;
????delta_q?=?theta（1）?*?（norminv（beta）?*?（2?*?theta（2）?*?h^2）^0.5?/?theta（1）?+?1?+?theta（2）?*?h?*?（h?-?1）?/?theta（1）^2）^（1/h）;
????
????xi?=?eye（n）;
????M_M?=?diag（sqrt_latent（1?:?i））?*?P‘?*?xi（:?1）;
????f_t（i）?=?2?*?（1?/?（svd（M_M）））?*?delta_t;
????
????C_C?=?tilde_P?*?tilde_P‘?*?xi（:?1）;
????f_q（i）?=?2?*?（1?/?（svd（C_C）））?*?delta_q;
end
f_t
[min_f_t?num_PC1]?=?min（f_t）
f_q
[min_f_q?num_PC2]?=?min（f_q）

Y?=?Creatmat（1000）;
W?=?zscore（Y）;
[m2?n2]?=?size（Y）;

%******加入故障******??????????????????????????????????
b?=?800;???????????????????????????????????????????????????????????????????%故障個(gè)數(shù)
%?h?=?randn（b?1）;?????????????????????????????????????????????????????????%加入干擾的類(lèi)型：精度下降干擾
h?=?linspace（12?12?b）‘;????????????????????????????????????????????????????%加入干擾的類(lèi)型：幅值為1的恒值干擾
W（m2?-?b?+?1?:?m2?1）?=?W（m2?-?b?+?1?:?m2?1）?+?h;?????????????????????????%疊加干擾
M?=?zscore（train_data）;
N?=?W;

%******建模******
S?=?cov（M）;????????????????????
[coeff?latent]?=?eig（S）;??????????
L?=?fliplr（diag（latent）‘）;?????
U?=?fliplr（coeff）;
sqrtL?=?sqrt（L）;
V?=?sqrtL（ones（size（M?2）?1）?:）?*?U;
%******PCA******
A1?=?num_PC1;???????????????????????????????????????????????????????????????????????????
P1?=?U（:?1?:?A1）;????????????????????????????????????????????????????????%取前A列作為新特征向量

Ts1?=?（diag（N?*?P1?*?diag（1./?L（1?:?A1））?*?P1‘?*?N‘））‘;????????????????????%計(jì)算T2統(tǒng)計(jì)量
I1?=?eye（size（P1?1））;
SPE1?=?（diag（N?*?（I1?-?P1?*?P1‘）?*?（I1?-?P1?*?P1‘）?*?N‘））‘;????????????????%計(jì)算SPE統(tǒng)計(jì)量

beta?=?0.99;
TsBeta1?=?A1?*?（m?-?1）?/?（m?-?A1）?*?finv（beta?A1?m?-?A1）;????????????????%求T2控制限??????????????????????????????????????????????????????????
theta1?=?zeros（3?1）;
for?i?=?1?:?3
????for?j?=?A1?+?1?:?size（M?2）
????????theta1（i）?=?theta1（i）?+?L（j）^（i）;
????end
end
h1?=?1?-?2?*?theta1（1）?*?theta1（3）?/?3?/?（theta1（2）^2）;
SPEBeta1?=?theta1（1）?*?（norminv（

?屬性????????????大小?????日期????時(shí)間???名稱(chēng)
-----------?---------??----------?-----??----
?????文件????????7643??2014-07-17?21:22??A_LinJie_PCA.m
?????文件?????????133??2014-01-29?18:02??Creatmat.m
?????文件?????????252??2014-01-29?18:02??creat_mat.m