1 x_avr=0;
2 y_avr=0;
3 x=1:1:50;
4 for i=1:1:50
5 y(i)=x(i)+2*(rand-1);
6 end
7 subplot(2,1,1);
8 plot(x,y,'*b');
9 axis([0 50 0 50 ]);
10 for i=1:50
11 x_avr=x_avr+x(i);
12 y_avr=y_avr+y(i);
13 end
14 x_avr=x_avr/50;
15 y_avr=y_avr/50;
16 for i=1:50
17 sig(1)=(x(i)-x_avr)*(x(i)-x_avr);
18 sig(2)=(x(i)-x_avr)*(y(i)-y_avr);
19 sig(3)=sig(2);
20 sig(4)=(y(i)-y_avr)*(y(i)-y_avr);
21 end
22 sig=sig/50;
23 sigma=[sigma(1),sig(2);sig(3),sig(4)];
24 d=[0 0 ; 0 0];
25 v=[0 0 ; 0 0];
26 [v,d]=eig(sigma);
27
28 %科学技术是第一生产力啊
29 %旋转数据
30 d
31 u2=[v(1,1);v(2,1)]%特征向量是列向量
32 u1=[v(1,2);v(2,2)]
33 for i=1:50
34 xr(i)=u1'*[x(i);y(i)];
35 yr(i)=u2'*[x(i);y(i)];
36 end
37 subplot(2,1,2);
38 plot(xr,yr,'*b');
39 axis([-100 100 -100 100])