统计建模与R软件课后习题答案25章.docx

上传人:b****3 文档编号:5450998 上传时间:2022-12-16 格式:DOCX 页数:19 大小:26.25KB
下载 相关 举报
统计建模与R软件课后习题答案25章.docx_第1页
第1页 / 共19页
统计建模与R软件课后习题答案25章.docx_第2页
第2页 / 共19页
统计建模与R软件课后习题答案25章.docx_第3页
第3页 / 共19页
统计建模与R软件课后习题答案25章.docx_第4页
第4页 / 共19页
统计建模与R软件课后习题答案25章.docx_第5页
第5页 / 共19页
点击查看更多>>
下载资源
资源描述

统计建模与R软件课后习题答案25章.docx

《统计建模与R软件课后习题答案25章.docx》由会员分享,可在线阅读,更多相关《统计建模与R软件课后习题答案25章.docx(19页珍藏版)》请在冰豆网上搜索。

统计建模与R软件课后习题答案25章.docx

统计建模与R软件课后习题答案25章

第二章答案:

Ex2.1

x<-c(1,2,3)

y<-c(4,5,6)

e<-c(1,1,1)

z=2*x+y+e

z1=crossprod(x,y)#z1为x1与x2的内积或者x%*%y

z2=tcrossprod(x,y)#z1为x1与x2的外积或者x%o%y

z;z1;z2

要点:

基本的列表赋值方法,内积和外积概念。

内积为标量,外积为矩阵。

Ex2.2

A<-matrix(1:

20,c(4,5));A

B<-matrix(1:

20,nrow=4,byrow=TRUE);B

C=A+B;C

#不存在AB这种写法

E=A*B;E

F<-A[1:

3,1:

3];F

H<-matrix(c(1,2,4,5),nrow=1);H

#H起过渡作用,不规则的数组下标

G<-B[,H];G

要点:

矩阵赋值方法。

默认是byrow=FALSE,数据按列放置。

取出部分数据的方法。

可以用数组作为数组的下标取出数组元素。

Ex2.3

x<-c(rep(1,times=5),rep(2,times=3),rep(3,times=4),rep(4,times=2));x#或者省略times=,如下面的形式

x<-c(rep(1,5),rep(2,3),rep(3,4),rep(4,2));x

要点:

rep()的使用方法。

rep(a,b)即将a重复b次

Ex2.4

n<-5;H<-array(0,dim=c(n,n))

for(iin1:

n){for(jin1:

n){H[i,j]<-1/(i+j-1)}};H

G<-solve(H);G#求H的逆矩阵

ev<-eigen(H);ev#求H的特征值和特征向量

要点:

数组初始化;for循环的使用

待解决:

如何将很长的命令(如for循环)用几行打出来再执行?

每次想换行的时候一按回车就执行了还没打完的命令...

Ex2.5

StudentData<-data.frame(name=c("zhangsan","lisi","wangwu","zhaoliu","dingyi"),sex=c("F","M","F","M","F"),age=c("14","15","16","14","15"),height=c("156","165","157","162","159"),weight=c("42","49","41.5","52","45.5"));StudentData

要点:

数据框的使用

待解决:

SSH登陆linux服务器中文显示乱码。

此处用英文代替。

Ex2.6

write.table(StudentData,file="studentdata.txt")

#把数据框StudentData在工作目录里输出,输出的文件名为studentdata.txt.

StudentData_a<-read.table("studentdata.txt");StudentData_a

#以数据框的形式读取文档studentdata.txt,存入数据框StudentData_a中。

write.csv(StudentData_a,"studentdata.csv")

#把数据框StudentData_a在工作目录里输出,输出的文件名为studentdata.csv,可用Excel打开.

要点:

读写文件。

read.table("file")

write.table(Rdata,"file")

read.csv("file")

write.csv(Rdata,"file")

外部文件,不论是待读入或是要写出的,命令中都得加双引号。

Ex2.7

Fun<-function(n){

if(n<=0)

list(fail="pleaseinputaintegerabove0!

")

else{

repeat{

if(n==1)break

elseif(n%%2==0){n<-n/2}

elsen<-3*n+1

}

list("sucess!

")

}

}

在linux下新建一个R文件,输入上述代码,保存为"2.7.R"

然后在当前目录下进入R环境,输入source("2.7.R"),即打开了这个程序脚本。

然后就可以执行函数了。

输入Fun(67),显示

"sucess!

"

输入Fun(-1),显示

$fail

[1]"pleaseinputaintegerabove0!

"

待解决:

source("*.R")是可以理解为载入这个R文件吧?

如何在R环境下关闭R文件呢?

Ex3.1

新建txt文件如下:

3.1.txt

74.379.575.073.575.874.073.567.275.873.578.875.673.575.075.8

72.079.576.573.579.568.875.078.872.068.876.573.572.775.070.4

78.078.874.364.376.574.374.770.472.776.570.472.075.875.870.4

76.565.077.273.572.780.572.065.080.371.277.676.568.873.577.2

80.572.074.369.781.267.381.667.372.784.369.774.371.274.375.0

72.075.467.381.675.071.271.269.773.570.475.072.767.370.376.5

73.572.068.073.568.074.372.772.774.370.4

编写一个函数(程序名为data_outline.R)描述样本的各种描述性统计量。

data_outline<-function(x){

n<-length(x)

m<-mean(x)

v<-var(x)

s<-sd(x)

me<-median(x)

cv<-100*s/m

css<-sum((x-m)^2)

uss<-sum(x^2)

R<-max(x)-min(x)

R1<-quantile(x,3/4)-quantile(x,1/4)

sm<-s/sqrt(n)

g1<-n/((n-1)*(n-2))*sum((x-m)^3)/s^3

g2<-((n*(n+1))/((n-1)*(n-2)*(n-3))*sum((x-m)^4)/s^4-(3*(n-1)^2)/((n-2)*(n-3)))

data.frame(N=n,Mean=m,Var=v,std_dev=s,Median=me,std_mean=sm,CV=cv,CSS=css,USS=uss,R=R,R1=R1,Skewness=g1,Kurtosis=g2,row.names=1)

}

进入R,

source("data_outline.R")#将程序调入内存

serumdata<-scan("3.1.txt");serumdata#将数据读入向量serumdata。

data_outline(serumdata)

结果如下:

NMeanVarstd_devMedianstd_meanCVCSSUSSR

110073.69615.416753.92641773.50.39264175.3278571526.258544636.320

R1SkewnessKurtosis

14.60.038542490.07051809

要点:

read.table()用于读表格形式的文件。

上述形式的数据由于第七行缺几个数据,故用read.table()不能读入。

scan()可以直接读纯文本文件。

scan()和matrix()连用还可以将数据存放成矩阵形式。

X<-matrix(scan("3.1.txt",0),ncol=10,byrow=TRUE)#将上述数据放置成10*10的矩阵。

scan()还可以从屏幕上直接输入数据。

Y<-scan()

然后按提示输入即可。

结束输入时按回车即可。

Ex3.2

>hist(serumdata,freq=FALSE,col="purple",border="red",density=3,angle=60,main=paste("thehistogramofserumdata"),xlab="age",ylab="frequency")#直方图。

col是填充颜色。

默认空白。

border是边框的颜色,默认前景色。

density是在图上画条纹阴影,默认不画。

angle是条纹阴影的倾斜角度(逆时针方向),默认45度。

main,xlab,ylab是标题,x和y坐标轴名称。

>lines(density(serumdata),col="blue")#密度估计曲线。

>x<-64:

85

>lines(x,dnorm(x,mean(serumdata),sd(serumdata)),col="green")#正态分布的概率密度曲线

>plot(ecdf(serumdata),verticals=TRUE,do.p=FALSE)#绘制经验分布图

>lines(x,pnorm(x,mean(serumdata),sd(serumdata)),col="blue")#正态经验分布

>qqnorm(serumdata,col="purple")#绘制QQ图

>qqline(serumdata,col="red")#绘制QQ直线

Ex3.3

>stem(serumdata,scale=1)#作茎叶图。

原始数据小数点后数值四舍五入。

Thedecimalpointisatthe|

64|300

66|23333

68|00888777

70|34444442222

72|0000000777777755555555555

74|033333333700000004688888

76|5555555226

78|0888555

80|355266

82|

84|3

>boxplot(serumdata,col="lightblue",notch=T)#作箱线图。

notch表示带有缺口。

>fivenum(serumdata)#五数总结

[1]64.371.273.575.884.3

Ex3.4

>shapiro.test(serumdata)#正态性Shapori-Wilk检验方法

Shapiro-Wilknormalitytest

data:

serumdata

W=0.9897,p-value=0.6437

结论:

p值>0.05,可认为来自正态分布的总体。

>ks.test(serumdata,"pnorm",mean(serumdata),sd(serumdata))#Kolmogrov-Smirnov检验,正态性

One-sampleKolmogorov-Smirnovtest

data:

serumdata

D=0.0701,p-value=0.7097

alternativehypothesis:

two-sided

Warningmessage:

Inks.test(serumdata,"pnorm",mean(serumdata),sd(serumdata)):

cannotcomputecorrectp-valueswithties

结论:

p值>0.05,可认为来自正态分布的总体。

注意,这里的警告信息,是因为数据中有重复的数值,ks检验要求待检数据时连续的,不允许重复值。

Ex3.5

>y<-c(2,4,3,2,4,7,7,2,2,5,4,5,6,8,5,10,7,12,12,6,6,7,11,6,6,7,9,5,5,10,6,3,10)#输入数据

>f<-factor(c(rep(1,11),rep(2,10),rep(3,12)))#因子分类

>plot(f,y,col="lightgreen")#plot()生成箱线图

>x<-c(2,4,3,2,4,7,7,2,2,5,4)

>y<-c(5,6,8,5,10,7,12,12,6,6)

>z<-c(7,11,6,6,7,9,5,5,10,6,3,10)

>boxplot(x,y,z,names=c("1","2","3"),col=c(5,6,7))#boxplot()生成箱线图

结论:

第2和第3组没有显著差异。

第1组合其他两组有显著差异。

Ex3.6

数据太多,懒得录入。

离散图应该用plot即可。

Ex3.7

>studata<-read.table("3.7.txt")#读入数据

>data.frame(studata)#转化为数据框

V1V2V3V4V5V6

11alicef1356.584.0

22beckaf1365.398.0

33gailf1464.390.0

44karenf1256.377.0

55kathyf1259.884.5

66maryf1566.5112.0

77sandyf1151.350.5

88sharonf1562.5112.5

99tammyf1462.8102.5

1010alfredm1469.0112.5

1111dukem1463.5102.5

1212guidom1567.0133.0

1313jamesm1257.383.0

1414jefferym1362.584.0

1515johnm1259.099.5

1616philipm1672.0150.0

1717robertm1264.8128.0

1818thomasm1157.585.0

1919williamm1566.5112.0

>names(studata)<-c("stuno","name","sex","age","height","weight"),studata#给各列命名

stunonamesexageheightweight

11alicef1356.584.0

22beckaf1365.398.0

33gailf1464.390.0

...

>attach(studata)#将数据框调入内存

>plot(weight~height,col="red")#体重对于身高的散点图

>coplot(weight~height|sex,col="blue")#不同性别,体重与身高的散点图

>coplot(weight~height|age,col="blue")#不同年龄,体重与身高的散点图

>coplot(weight~height|age+sex,col="blue")#不同年龄和性别,体重与身高的散点图

Ex3.8

>x<-seq(-2,3,0.05)

>y<-seq(-1,7,0.05)

>f<-function(x,y)x^4-2*x^2*y+x^2-2*x*y+2*y^2+4.5*x-4*y+4

>z<-outer(x,y,f)#必须做外积运算才能绘出三维图形

>contour(x,y,z,levels=c(0,1,2,3,4,5,10,15,20,30,40,50,60,80,100),col="blue")#二维等值线

>persp(x,y,z,theta=120,phi=0,expand=0.7,col="lightblue")#三位网格曲面

Ex3.9

>attach(studata)

>cor.test(height,weight)#Pearson相关性检验

Pearson'sproduct-momentcorrelation

data:

heightandweight

t=7.5549,df=17,p-value=7.887e-07

alternativehypothesis:

truecorrelationisnotequalto0

95percentconfidenceinterval:

0.70443140.9523101

sampleestimates:

cor

0.8777852

由此可见身高和体重是相关的。

Ex4.2

指数分布,λ的极大似然估计是n/sum(Xi)

>x<-c(rep(5,365),rep(15,245),rep(25,150),rep(35,100),rep(45,70),rep(55,45),rep(65,25))

>lamda<-length(x)/sum(x);lamda

[1]0.05

Ex4.3

Poisson分布P(x=k)=λ^k/k!

*e^(-λ)

其均数和方差相等,均为λ,其含义为平均每升水中大肠杆菌个数。

取均值即可。

>x<-c(rep(0,17),rep(1,20),rep(2,10),rep(3,2),rep(4,1))

>mean(x)

[1]1

平均为1个。

Ex4.4

>obj<-function(x){f<-c(-13+x[1]+((5-x[2])*x[2]-2)*x[2],-29+x[1]+((x[2]+1)*x[2]-14)*x[2]);sum(f^2)}#其实我也不知道这是在干什么。

所谓的无约束优化问题。

>x0<-c(0.5,-2)

>nlm(obj,x0)

$minimum

[1]48.98425

$estimate

[1]11.4127791-0.8968052

$gradient

[1]1.411401e-08-1.493206e-07

$code

[1]1

$iterations

[1]16

Ex4.5

>x<-c(54,67,68,78,70,66,67,70,65,69)

>t.test(x)#t.test()做单样本正态分布区间估计

OneSamplet-test

data:

x

t=35.947,df=9,p-value=4.938e-11

alternativehypothesis:

truemeanisnotequalto0

95percentconfidenceinterval:

63.158571.6415

sampleestimates:

meanofx

67.4

平均脉搏点估计为67.4,95%区间估计为63.158571.6415。

>t.test(x,alternative="less",mu=72)#t.test()做单样本正态分布单侧区间估计

OneSamplet-test

data:

x

t=-2.4534,df=9,p-value=0.01828

alternativehypothesis:

truemeanislessthan72

95percentconfidenceinterval:

-Inf70.83705

sampleestimates:

meanofx

67.4

p值小于0.05,拒绝原假设,平均脉搏低于常人。

要点:

t.test()函数的用法。

本例为单样本;可做双边和单侧检验。

Ex4.6

>x<-c(140,137,136,140,145,148,140,135,144,141);x

[1]140137136140145148140135144141

>y<-c(135,118,115,140,128,131,130,115,131,125);y

[1]135118115140128131130115131125

>t.test(x,y,var.equal=TRUE)

TwoSamplet-test

data:

xandy

t=4.6287,df=18,p-value=0.0002087

alternativehypothesis:

truedifferenceinmeansisnotequalto0

95percentconfidenceinterval:

7.5362620.06374

sampleestimates:

meanofxmeanofy

140.6126.8

期望差的95%置信区间为7.5362620.06374。

要点:

t.test()可做两正态样本均值差估计。

此例认为两样本方差相等。

ps:

我怎么觉得这题应该用配对t检验?

Ex4.7

>x<-c(0.143,0.142,0.143,0.137)

>y<-c(0.140,0.142,0.136,0.138,0.140)

>t.test(x,y,var.equal=TRUE)

TwoSamplet-test

data:

xandy

t=1.198,df=7,p-value=0.2699

alternativehypothesis:

truedifferenceinmeansisnotequalto0

95percentconfidenceinterval:

-0.0019963510.006096351

sampleestimates:

meanofxmeanofy

0.141250.13920

期望差的95%的区间估计为-0.0019963510.006096351

Ex4.8

接Ex4.6

>var.test(x,y)

Ftesttocomparetwovariances

data:

xandy

F=0.2353,numdf=9,denomdf=9,p-value=0.04229

alternativehypothesis:

trueratioofvariancesisnotequalto1

95percentconfidenceinterval:

0.058452760.94743902

sampleestimates:

ratioofvariances

0.2353305

要点:

var.test可做两样本方差比的估计。

基于此结果可认为方差不等。

因此,在Ex4.6中,计算期望差时应该采取方差不等的参数。

>t.test(x

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 医药卫生 > 基础医学

copyright@ 2008-2022 冰豆网网站版权所有

经营许可证编号:鄂ICP备2022015515号-1