附录
Clear
clc
file=dir('C:\Users\孙晓慧\Documents\MATLAB\normal*.txt');
for i=1:length(file)
temp=load(['C:\Users\孙晓慧\Documents\MATLAB\normal',file(i).name]);
eval([file(i).name(1:end-4),'=temp;'])
normal=temp;
normal(1,:)=[];
normal(:,18)=[];
normal(:,17)=[];
normal(:,16)=[];
normal(:,15)=[];
normal(:,8)=[];
normal(:,5)=[];
normal(:,4)=[];
normal(:,1)=[];
[m,n]=size(normal);
twodata=[];for p=1:399
for j=1:n
twodata=[twodata normal(p,j)];
end
end
pensim_saw_data(i,:)= twodata;
end
normaldata=pensim_saw_data(1:50,:);
%随机提取12批正常数据作为校验数据
validation=[normaldata(2,:);normaldata(4,:);normaldata(6,:);normaldata(8,:);
normaldata(12,:);normaldata(16,:);normaldata(20,:);normaldata(25,:);
normaldata(28,:);normaldata(30,:);normaldata(35,:);normaldata(48,:)];
%提取校验数据后，在正常数据矩阵中删掉
normaldata(48,:)=[];normaldata(35,:)=[];normaldata(30,:)=[];normaldata(28,:)=[];
normaldata(25,:)=[];normaldata(20,:)=[];normaldata(16,:)=[];normaldata(12,:)=[];
normaldata(8,:)=[];normaldata(6,:)=[];normaldata(4,:)=[];normaldata(2,:)=[];
file=dir('C:\Users\孙晓慧\Documents\MATLAB\fault*.txt');
for i=1:length(file)
temp=load(['C:\Users\孙晓慧\Documents\MATLAB\fault',file(i).name]);
eval([file(i).name(1:end-4),'=temp;'])
normal=temp;
normal(1,:)=[];
normal(:,18)=[];
normal(:,17)=[];
normal(:,16)=[];
normal(:,15)=[];
normal(:,8)=[];
normal(:,5)=[];
normal(:,4)=[];
normal(:,1)=[];
[m,n]=size(normal);
twodata=[];
for p=1:399
for j=1:n
twodata=[twodata normal(p,j)];
end
end
pensim_saw_data1(i,:)= twodata;
end
pt=pensim_saw_data1;
pt_test=[validation;pt];%112行为正常校验数据，1330行为故障数据
[normalindex1,mnormal,stsnormal]=auto(normaldata);%建模矩阵标准化
[neighbors distances] = kNearestNeighbors(normalindex1,normalindex1,38);
C=distances;
for i=1:38
DIS(i,1)=C(i,2)+C(i,3)+C(i,4)+C(i,5)%+C(i,6);
end
[Df,x]=ksdensity(DIS(:,1),'function','cdf');
plot(x,Df)
for k=1:100
Dfmin95=min(abs(Df-0.95));
if(abs(Df(1,k)-0.95)==Dfmin95)
d=k;
contr95=x(1,d);
break
end
end
d1=1;
for k=1:100
Dfmin99=min(abs(Df-0.99));
if(abs(Df(1,k)-0.99)==Dfmin99)
d1=k;
contr99=x(1,d1);
break
end
end
M=pt_test; %包含有校验数据和故障数据的数据文件，32（8517）
a=mnormal;
c=stsnormal;
[p,q]=size(M);
for j=1:q
for i=1:p
M(i,j)=(M(i,j)-a(1,j))/c(1,j);
end
end
[neighbors1 distances1] = kNearestNeighbors(normalindex1,M,34);
C1=distances1;
for i=1:34
DIS2(i,1)=C1(i,1)+C1(i,2)+C1(i,3)+C1(i,4)%+C1(i,5);
end
Validation=DIS2(1:12,:);
fault=DIS2(13:34,:);
figure
plot(contr95ones(1,80),'k--')
hold on
i=1:38
plot(i,DIS,'ko');
hold on;
j=1:12
plot(38+j,Validation,'bs');
hold on
for iii=1:22
plot(50+iii,fault(iii,1),'r');text(50+iii,fault(iii,1),num2str(iii));
end
hold on
%legend('Threshold','Training','Validation-normal','Validation-fault',2)
xlabel('batch')
ylabel('knn squared distance(k=3)')
title(' pen_knn')
hold off
function [ax,mx,stdx] = auto(x)
%AUTO Autoscales matrix to mean zero unit variance
% Autoscales a matrix (x) and returns the resulting matrix (ax)
% with mean-zero unit variance columns, a vector of means (mx)
% and a vector of standard deviations (stdx) used in the scaling.
%
%I/O: [ax,mx,stdx] = auto(x);
%
%See also: MDAUTO, MDMNCN, MDRESCAL, MDSCALE, MNCN, SCALE, RESCALE
%Copyright Eigenvector Research, Inc. 1991-98
%Modified 11/93
%Checked on MATLAB 5 by BMW 1/4/97
[m,n] = size(x);
mx = mean(x);
stdx = std(x);
ax = (x-mx(ones(m,1),:))./stdx(ones(m,1),:);
End
function [neighborIds neighborDistances] = kNearestNeighbors(dataMatrix, queryMatrix, k)
%--------------------------------------------------------------------------
% Program to find the k - nearest neighbors (kNN) within a set of points.
% Distance metric used: Euclidean distance
%
% Usage:
% [neighbors distances] = kNearestNeighbors(dataMatrix, queryMatrix, k);
% dataMatrix (N x D) - N vectors with dimensionality D (within which we search for the nearest neighbors)
% queryMatrix (M x D) - M query vectors with dimensionality D
% k (1 x 1) - Number of nearest neighbors desired
%
% Example:
% a = [1 1; 2 2; 3 2; 4 4; 5 6];
% b = [1 1; 2 1; 6 2];
% [neighbors distances] = kNearestNeighbors(a,b,2);
%
% Output:
% neighbors =
% 1 2
% 1 2
% 4 3
%
% distances =
% 0 1.4142
% 1.0000 1.0000
% 2.8284 3.0000
%--------------------------------------------------------------------------
neighborIds = zeros(size(queryMatrix,1),k);
neighborDistances = neighborIds;
numDataVectors = size(dataMatrix,1);
numQueryVectors = size(queryMatrix,1);
for i=1:numQueryVectors,
dist = sum((repmat(queryMatrix(i,:),numDataVectors,1)-dataMatrix).^2,2);
[sortval sortpos] = sort(dist,'ascend');
neighborIds(i,:) = sortpos(1:k);
% neighborDistances(i,:) = sqrt(sortval(1:k));
neighborDistances(i,:) = sortval(1:k);
end
end
