南京师范大学泰州学院本科毕业论文

附录

附录一：取样点分布图

A = xlsread('cumcm2011A附件_数据.xls','排序后取样点表','B4:D322'); B = xlsread('cumcm2011A附件_数据.xls','排序后浓度表','B4:I322'); x=A(:,1);y=A(:,2);z=B(:,7); x1=x(1:44);y1=y(1:44);z1=z(1:44); plot3(x1,y1,z1,'*b') hold on

x2=x(45:80);y2=y(45:80);z2=z(45:80); scatter3(x2,y2,z2,'pg')

x3=x(81:146);y3=y(81:146);z3=z(81:146); scatter3(x3,y3,z3,'+k')

x4=x(147:284);y4=y(147:284);z4=z(147:284); scatter3(x4,y4,z4,'.r')

x5=x(285:319);y5=y(285:319);z5=z(285:319); scatter3(x5,y5,z5,'hm') xlabel('x轴 单位（m）') ylabel('y轴 单位（m）') zlabel('海拔 单位(m）') title('取样点分布图')

legend('生活区取样点','工业区取样点','山区取样点','主干道路区取样点','公园绿地区取样点') x=A(:,1);y=A(:,2);z=A(:,3);

[X,Y,Z]=griddata(x,y,z,linspace(min(x),max(x))',linspace(min(y),max(y)),'v4');%插值 meshc(X,Y,Z)

附录二：取样点等高分布图：

A = xlsread('cumcm2011A附件_数据.xls','排序后取样点表','B4:D322'); x=A(:,1);y=A(:,2);z=A(:,3);

[X,Y,Z]=griddata(x,y,z,linspace(min(x),max(x))',linspace(min(y),max(y)),'v4');%插值 contour(X,Y,Z,30) hold on

x1=x(1:44);y1=y(1:44);z1=z(1:44); plot(x1,y1,'*b')

x2=x(45:80);y2=y(45:80);z2=z(45:80); plot(x2,y2,'pg')

x3=x(81:146);y3=y(81:146);z3=z(81:146); plot(x3,y3,'+k')

x4=x(147:284);y4=y(147:284);z4=z(147:284); plot(x4,y4,'.r')

x5=x(285:319);y5=y(285:319);z5=z(285:319); plot(x5,y5,'hm')

xlabel('x轴 单位（m）') ylabel('y轴 单位（m）') title('取样点等高分布图')

32

南京师范大学泰州学院本科毕业论文

legend('海拔等高线','生活区取样点','工业区取样点','山区取样点','主干道路区取样点','公园绿地区取样点') x=A(:,1);y=A(:,2);z=A(:,3);

[X,Y,Z]=griddata(x,y,z,linspace(min(x),max(x))',linspace(min(y),max(y)),'v4');%插值 contour(X,Y,Z,30)

附录三：工业区As、Cd、Cr、Cu元素的浓度与海拔关系图 A = xlsread('cumcm2011A附件_数据.xls','工业区','A2：A37'); B = xlsread('cumcm2011A附件_数据.xls','工业区','B2:I37'); subplot(2,2,1) As = B(:,1); plot(A,As)

xlabel('海拔 单位（m）') ylabel('As浓度(μg/g)')

title('工业区As元素浓度与海拔关系折线图') subplot(2,2,2) Cd = B(:,2); plot(A,Cd)

xlabel('海拔 单位（m）') ylabel('Cd浓度(ng/g)')

title('工业区Cd元素浓度与海拔关系折线图') subplot(2,2,3) Cr = B(:,3); plot(A,Cr)

xlabel('海拔 单位（m）') ylabel('Cr浓度(μg/g)')

title('工业区Cr元素浓度与海拔关系折线图') subplot(2,2,4) Cu = B(:,4); plot(A,Cu)

xlabel('海拔 单位（m）') ylabel('Cu浓度(μg/g)')

title('工业区Cu元素浓度与海拔关系折线图')

工业区Hg、Ni、Pb、Zn元素的浓度与海拔关系图，山区As、Cd、Cr、Cu元素的浓度与海拔关系，山区Hg、Ni、Pb、Zn元素的浓度与海拔关系图程序类似略 附录四：因子分析程序

x = xlsread('cumcm2011A附件_数据.xls','附件2','B4:I322'); x = zscore(x)%数据标准化 r = cov(x)%相关矩阵 [vec,val,con] = pcacov(r) val,con

f1 = repmat(sign(sum(vec)),size(vec,1),1); vec = vec.*f1;

f2 = repmat(sqrt(val)',size(vec,1),1); a = vec.*f2%载荷矩阵

[b,t] = rotatefactors(a(:,1:3),'method','varimax')

33

南京师范大学泰州学院本科毕业论文

bz = [b,a(:,4:end)] gx = sum(bz.^2) gxv = gx/sum(gx) 附录五：变量聚类分析 a = [1.0000 0 0.2547 0.1890 0.1597 0.0644 0.3166 0.2899 0.2469 for i = 1:8 a(i,i) = 0; end b = a(:); b = nonzeros(b); b = b'; b = 1-b;

z = linkage(b,'complete');%最大距离 y = cluster(z,2)%阈值为2 dendrogram(z)%绘图 % ind1 = find(y==2); % ind1 = ind1' % ind2 = find(y==1); % ind2 = ind2'

title('土壤重金属相关元素聚类分析图')

set(gca,'XTicklabel',['Cr'; 'Ni';'Cu';'Cd';'Pb';'Zn','As';'Hg']) % h = get(gca)

% set(h,'XTickLabel',['Cr';'Ni';'Cu';'Cd';'Pb';'Zn';'As';'Hg']) % set(gca,'XTickLabel',['Cr Ni Cu Cd Pb Zn As Hg']) 附录六：熵值法计算权重 x = [6.27 289.9613636

237.0086364;

53.40916667 38.95969697 58.05391304 43.636

127.5358333 17.31727273 62.21492754

642.3552778 40.95606061 446.8225362

19.81166667 15.45378788 17.61710145 15.28971429

93.04083333 36.55590909 63.5342029 60.70857143

277.9275 73.29424242 242.8548551 154.2422857];

7.25 393.1111111 4.04 152.319697 5.71 360.0144928 6.26 280.5428571 [n,m] = size(x) % x =load('sj.txt') %标准化处理 x = abs(zscore(x)) %计算比重 s=sum(x)

69.01840909

49.40318182

93.04068182

18.34227273

69.10636364

0

0

0 0

0 0

0 0 0

0 0 0

0 0 0

0 0 0

0 0 0

0 0

0 0

0

0 1.0000];

1.0000 0.4937

1.0000 0.3524 0.3967 0.2647 0.3294 0.6603 0.4312

1.0000 0.5316 0.1032 0.7158 0.3828 0.4243

1.0000 0.4167 0.4946 0.5200 0.3873

1.0000 0.1029 0.2981 0.1958

1.0000 0.3068 0.4364

30.19171429 114.9917143

34

南京师范大学泰州学院本科毕业论文

p1 = x(1,:)./s p2 = x(2,:)./s p3 = x(3,:)./s p4 = x(4,:)./s p5 = x(5,:)./s p = [p1;p2;p3;p4;p5] %计算熵值 k = log(1/n)

e1= k*(p1(:,1)+p2(:,1)+p3(:,1)+p4(:,1)+p5(:,1)+log(p1(:,1))+log(p2(:,1))+log(p3(:,1))+log(p4(:,1))+log(p4(:,1))) e2= k*(p1(:,2)+p2(:,2)+p3(:,2)+p4(:,2)+p5(:,2)+log(p1(:,2))+log(p2(:,2))+log(p3(:,2))+log(p4(:,2))+log(p4(:,2))) e3= k*(p1(:,3)+p2(:,3)+p3(:,3)+p4(:,3)+p5(:,3)+log(p1(:,3))+log(p2(:,3))+log(p3(:,3))+log(p4(:,3))+log(p4(:,3))) e4= k*(p1(:,4)+p2(:,4)+p3(:,4)+p4(:,4)+p5(:,4)+log(p1(:,4))+log(p2(:,4))+log(p3(:,4))+log(p4(:,4))+log(p4(:,4))) e5= k*(p1(:,5)+p2(:,5)+p3(:,5)+p4(:,5)+p5(:,5)+log(p1(:,5))+log(p2(:,5))+log(p3(:,5))+log(p4(:,5))+log(p4(:,5))) e6= k*(p1(:,6)+p2(:,6)+p3(:,6)+p4(:,6)+p5(:,6)+log(p1(:,6))+log(p2(:,6))+log(p3(:,6))+log(p4(:,6))+log(p4(:,6))) e7= k*(p1(:,7)+p2(:,7)+p3(:,7)+p4(:,7)+p5(:,7)+log(p1(:,7))+log(p2(:,7))+log(p3(:,7))+log(p4(:,7))+log(p4(:,7))) e8= k*(p1(:,8)+p2(:,8)+p3(:,8)+p4(:,8)+p5(:,8)+log(p1(:,8))+log(p2(:,8))+log(p3(:,8))+log(p4(:,8))+log(p4(:,8))) %计算差异系数

E = e1+e2+e3+e4+e5+e6+e7+e8; g1 = (1-e1)/(m-E) g2 = (1-e2)/(m-E) g3 = (1-e3)/(m-E) g4 = (1-e4)/(m-E) g5 = (1-e5)/(m-E) g6 = (1-e6)/(m-E) g7 = (1-e7)/(m-E) g8 = (1-e8)/(m-E)

% g1+g2+g3+g4+g5+g6+g7+g8 %计算权重

w = [g1,g2,g3,g4,g5,g6,g7,g8] %计算变量分值

xx = [3 3 3 4 3 2 3 4;3 4 3 5 5 2 4 4;2 2 2 2 2 2 2 2;3 3 3 5 5 2 3 4;3 3 2 3 3 2 3 3] w = 10*[w;w;w;w;w] V = xx.*w VV = sum(V,2) N = 10./VV

35