《数值计算方法》课程设计报告.docx
《《数值计算方法》课程设计报告.docx》由会员分享,可在线阅读,更多相关《《数值计算方法》课程设计报告.docx(21页珍藏版)》请在冰豆网上搜索。
《数值计算方法》课程设计报告
《数据分析方法》
课程实验报告
项目名称:
回归分析
学生姓名:
滕伟
学生学号:
200707020118
指导教师:
白林
完成日期:
2009年11月13日
1.实验内容
(1)掌握回归分析的思想和计算步骤;
(2)编写程序完成回归分析的计算,包括后续的显著性检验、残差分析、Box-Cox变换等内容。
2.模型建立与求解(数据结构与算法描述)
3.实验数据与实验结果
解:
根据所建立的模型在MATLAB中输入程序(程序见附录)得到以下结果:
(1)回归方程为:
说明该化妆品的消量和该城市人群收入情况关系不大,轻微影响,与使用该化妆品的人数有关。
的无偏估计:
(2)方差分析表如下表:
方差来源
自由度
平方和
均方
值
回归(
)
2
53845
26922
56795
2.28
误差(
)
12
56.883
4.703
总和(
)
14
53902
从分析表中可以看出:
值远大于
的值。
所以回归关系显著。
复相关
,所以回归效果显著。
解:
根据所建立的模型,在MATLAB中输入程序(程序见附录)得到如下结果:
(1)回归方程为:
在MTLAB中计算学生化残差(见程序清单二),所得到的学生化残差r的值由残差可知得到的r的值在(-1,1)的概率为0.645,在(-1.5,1.5)的概率为0.871,在(-2,2)之间的概率为0.968.
而服从正态分布的随机变量取值在(-1,1)之间的概率为0.68,在(-1.5,1.5)之间的概率为0.87,在(-2.2)之间的概率为0.95,所以相差较大,所以残差分析不合理,需要对数据变换。
取
=0.6进行Box-Cox变换
在MATLAB中输入程序(见程序代码清单二)
取
,所以得到r的值(r的值见附录二)其值在(-1,1)之间的个数大约为20/31=0.65,大致符合正态分布,所以重新拟合为:
拟合函数为:
通过F值,R值可以检验到,回归效果显著
(3)某医院为了了解病人对医院工作的满意程度
和病人的年龄
,病情的严重程度
和病人的忧虑程度
之间的关系,随机调查了该医院的23位病人,得数据如下表:
年龄
病情程度
忧虑程度
满意程度
50
51
2.3
48
36
46
2.3
57
40
48
2.2
66
41
44
1.8
70
28
43
1.8
89
49
54
2.9
36
42
50
2.2
46
45
48
2.4
54
52
62
2.9
26
29
50
2.1
77
29
48
2.4
89
43
53
2.4
67
38
55
2.2
47
34
51
2.3
51
53
54
2.2
57
36
49
2.0
66
33
56
2.5
79
29
46
1.9
88
33
49
2.1
60
55
51
2.4
49
29
52
2.3
77
44
58
2.9
52
43
50
2.3
60
(1)拟合线性回归模型
,通过残差分析与考察模型及有关误差分布正态性假定的合理性;
(2)若
(1)中模型合理,分别在
,
,
准则下选择最优回归方程,各准则下的选择结果是否一致?
(3)对
,用逐步回归法选择最优回归方程,其结果和
(2)中的数否一致?
(4)对选择的最优回归方程作残差分析,与
(1)中的相应结果比较,有何变化?
习题2.6
解:
(1)回归参数的
的最小二乘估计为:
。
在MATLAB中输入程序(见程序代码清单二)可得:
,
,
所以回归方程为:
对数据做Box-Cox变换,(由于
的取值在能力范围不好确定,所以经测试,取
=0.6进行Box-Cox变换
在MATLAB中输入程序(见程序代码清单二)
取
,所以得到r的值(r的值见附录二)其值在(-1,1)之间的个数大约为20/31=0.65,大致符合正态分布,所以重新拟合为:
拟合函数为:
通过F值,R值可以检验到,回归效果显著
习题2.9
解:
根据所建立的模型,在MATLAB中输入程序,得到以下结果:
(1)所得到的回归方程为:
(2)所得到的学生化残差见附录,通过对残差的分析,很明显不符合正态分布所以
(1)中所建立的模型不合理。
4.程序代码清单:
习题2.4
x=[12742450
11803254
13753802
12052838
1862347
12653782
1983008
13302450
11952137
1532560
14304020
13724427
12362660
11572088
13702605];
y=[162
120
223
131
67
169
81
192
116
55
252
232
144
103
212];
n=15;p=3
b=inv(x'*x)*x'*y
h=x*inv(x'*x)*x';
sse=y'*(eye(n,n)-h)*y
d2=1/(n-p)*y'*(eye(n,n)-h)*y
sst=y'*(eye(n,n)-(1/n)*ones(n,n))*y
ssr=y'*(h-1/n*ones(n,n))*y
msr=ssr/(p-1)
mse=sse/(n-p)
f=msr/mse
r2=1-sse/sst
习题2.6
x=[18.370
18.665
18.863
110.572
110.781
110.883
111.066
111.075
111.180
111.275
111.379
111.476
111.476
111.769
112.075
112.974
112.985
113.386
113.771
113.864
114.078
114.280
114.574
116.072
116.377
117.381
117.582
117.980
118.080
118.080
120.687];
y=[10.3
10.3
10.2
16.4
18.8
19.7
15.6
18.2
22.6
19.9
24.2
21.0
21.4
21.3
19.1
22.2
33.8
27.4
25.7
24.9
34.5
31.7
36.3
38.3
42.6
55.4
55.7
58.3
51.5
51.0
77.0];
n=31;p=3;
b=inv(x'*x)*x'*y;
h=x*inv(x'*x)*x';
sst=y'*(eye(n,n)-(1/n)*ones(n,n))*y
sse=y'*(eye(n,n)-h)*y
mse=sse/(n-p)
ssr=y'*(h-1/n*ones(n,n))*y
msr=ssr/(p-1)
f=msr/mse
r2=1-sse/sst
fori=1:
n
a=h(2*(i-1)+i)
end
t=sqrt((mse-mse*a))
q=y-(-57.9877+4.7082*x(:
2:
2)+0.3393*x(:
3:
3))
r=q/t
程序三
x=[18.370
18.665
18.863
110.572
110.781
110.883
111.066
111.075
111.180
111.275
111.379
111.476
111.476
111.769
112.075
112.974
112.985
113.386
113.771
113.864
114.078
114.280
114.574
116.072
116.377
117.381
117.582
117.980
118.080
118.080
120.687];
y=[10.3
10.3
10.2
16.4
18.8
19.7
15.6
18.2
22.6
19.9
24.2
21.0
21.4
21.3
19.1
22.2
33.8
27.4
25.7
24.9
34.5
31.7
36.3
38.3
42.6
55.4
55.7
58.3
51.5
51.0
77.0];
n=31;p=3;
m=0:
0.01:
1;
y=(y.^m-1)/m
b=inv(x'*x)*x'*y
h=x*inv(x'*x)*x';
sse=y'*(eye(n,n)-h)*y
mse=sse/(n-p)
f=msr/mse
r2=1-sse/sst
fori=1:
n
a=h(2*(i-1)+i)
end
t=sqrt((mse-mse*a))
q=y-(-57.9877+4.7082*x(:
2:
2)+0.3393*x(:
3:
3))
r=q/t
习题2.9
a=[150512.348
136462.357
140482.266
141441.870
128431.889
149542.936
142502.246
145482.454
152622.926
129502.177
129482.489
143532.467
138552.247
134512.351
153542.257
136492.066
133562.579
129461.988
133492.160
155512.449
129522.377
144582.952
143502.360]
y=a(:
5:
5)
x=a(:
1:
4)
n=23;p=4;
b=inv(x'*x)*x'*y
h=x*inv(x'*x)*x';
sst=y'*(eye(n,n)-(1/n)*ones(n,n))*y
sse=y'*(eye(n,n)-h)*y
mse=sse/(n-p)
ssr=y'*(h-1/n*ones(n,n))*y
msr=ssr/(p-1)
f=msr/mse
r2=1-sse/sst
fori=1:
n
a=h(2*(i-1)+i)
end
t=sqrt((mse-mse*a))
q=y-(162.8575-1.2103*x(:
2:
2)-0.6659*x(:
3:
3)-8.613*x(:
4:
4))
r=q/t
附录:
习题2.6学生化残差
r=
1.3857
1.4578
1.3656
0.1325
-0.2725
-0.3358
-0.1514
-0.2665
0.3002
-0.0740
0.5535
-0.1200
-0.0184
0.2006
-1.2333
-1.4358
0.5614
-1.6275
-1.2451
-0.9648
0.0273
-1.0948
0.2312
-0.8816
-0.5793
1.1303
0.8813
1.2355
-0.6102
-0.7372
2.1526
习题2.9学生化残差
r=
-0.0558
-1.1563
0.2408
0.1530
0.4069
-0.6459
-1.3453
-0.1718
-0.7476
0.0609
1.3545
1.1860
-1.3953
-1.6533
1.2882
-0.3350
1.4551
0.7065
-1.1911
0.7166
0.3590
0.5853
0.2236
上课纪律(20%)
实验过程及结果(40%)
实验报告质量(40%)
总分:
教师签字:
1.实验内容
(1)掌握主成份分析与典型相关分析的思想和计算步骤;
(2)编写程序完成主成份分析与典型相关性分析的计算;
2.模型建立与求解(数据结构与算法描述)
1.计算样本主成分的步骤:
(1)计算样本协方差矩阵S和相关系数矩阵R:
(2)计算S的特征值和相应的正交化特征向量:
,
(3)第K个样本的得分样本方差:
(4)前M个样本主成分的累加贡献率:
(5)选取m(m
2.计算样本典型变量相关系数的步骤:
(1)计算样本的协方差矩阵:
(2)计算A,B矩阵的特征值和正交化向量
(3)第K个样本典型相关变量为:
3.实验数据与实验结果
习题4.5
解:
在MATLAB中输入程序(见附录)
样本相关系数矩阵R为:
1
0.3336
-0.0545
-0.0613
-0.2894
0.1988
0.3487
0.3187
0.3336
1
-0.0229
0.3989
-0.1563
0.7111
0.4136
0.835
-0.0545
-0.0229
1
0.5333
0.4968
0.0328
-0.1391
-0.2584
-0.0613
0.3989
0.5333
1
0.6984
0.4679
-0.1713
0.3128
-0.2894
-0.1563
0.4968
0.6984
1
0.2801
-0.2083
-0.0812
0.1988
0.7111
0.0328
0.4679
0.2801
1
0.4168
0.7016
0.3487
0.4136
-0.1391
-0.1713
-0.2083
0.4168
1
0.3989
0.3187
0.835
-0.2584
0.3128
-0.0812
0.7016
0.3989
1
对应的特征值为:
3.0963
2.3672
0.92
0.7059
0.4984
0.0515
0.1308
0.2299
所以各主成分的贡献率为:
X1
0.387
X5
0.0623
X2
0.2959
X6
0.0064
X3
0.115
X7
0.0163
X4
0.0882
X8
0.0287
前两个主成分的累加贡献率为:
0.3870+0.2959=0.6859
各省市按照第一主成分排序,结果如下:
海南
河南
宁夏
西藏
广西
广东
陕西
湖北
辽宁
江苏
天津
内蒙古
山西
北京
四川
福建
甘肃
上海
黑龙江
新疆
青海
河北
吉林
浙江
湖南
云南
山东
安徽
贵州
江西
习题4.10
解:
在MATLAB中输入程序(程序见清单二):
得到相关系数矩阵R:
1
0.9362
0.4934
0.9362
1
0.7677
0.4934
0.7677
1
对应的特征值为:
0
0.4166
0.9091
4程序清单:
清单一
a=[8.3523.537.518.6217.42101.0411.21
9.2523.756.619.1917.7710.481.7210.51
8.1930.54.729.7816.287.62.5210.32
7.7329.25.429.4319.298.492.5210
9.4227.938.28.1416.179.421.559.76
9.1627.989.019.3215.999.11.8211.35
10.0628.6410.5210.0516.188.391.9610.81
9.0928.127.49.6217.2611.122.4912.65
9.4128.25.7710.816.3611.561.5312.17
8.728.127.2110.5319.4513.31.6611.96
6.9329.854.549.4916.6210.651.8813.61
8.6736.057.317.7516.6711.682.3812.88
9.9837.697.018.9416.1511.080.8311.67
6.7738.696.018.8214.7911.441.7413.23
8.1437.759.618.4913.159.761.2811.28
7.6735.718.048.3115.137.761.4113.25
7.939.778.4912.9419.2711.052.0413.29
7.1840.917.328.9417.612.751.1414.8
8.8233.77.5910.9818.8214.731.7810.1
6.2535.024.726.2810.037.151.9310.39
10.652.417.79.9812.5311.72.3114.69
7.2752.653.849.1613.0315.261.9814.57
13.4555.855.57.459.559.522.2116.3
10.8544.687.3214.5117.1312.081.2611.57
7.2145.797.6610.3616.5612.862.2511.69
7.6850.3711.3513.319.2514.592.7514.87
7.7848.44820.5122.1215.731.1516.61
7.9439.6520.9720.8222.5212.411.757.9
8.2864.34822.2220.0615.120.7222.89
12.4776.395.5211.2414.52225.4625.5];
r=corrcoef(a);
b=eig(r)
fori=1:
8
e=b(i)/sum(b)
end
清单二:
a=[606962976998
56538410378107
8069766699130
5580908085114
62756811613091
746470109101103
64716677102130
737064115110109
6867757685119
69827472133127
606761130134121
707478150158100
667478150131142
8370749998105
68669011985109
78637516498138
1037777160117121
77687414471153
667768778289
70707211493122
7565717770109
917493118115150
667573170147121
758276153132115
747166143105100
767064114113129
74908673106116
7477801168177
677169638770
78758010513280
6466718394133
718076818786
6375731208959
9010374107109101
6076619911198
48777511312497
669397136112122
74707610988105
607471729071
63756613010190
668086130117144
7767748392107
7067100150142146
737681119120119
789077122155149
73688010290122
7283681046996
6560701199489
5270769294100];
b=a';
r=corrcoef(b);
r11=r(1:
3,1:
3);
r21=r(4:
6,1:
3);
r12=r21';
r22=r(4:
6,4:
6);
R=corrcoef(inv(r11)*r12*inv(r22)*r21)
lamda=eig(R);
p=sqrt(lamda)