R语言的流行病学应用 课件 第6章 基本统计分析方法

第6章

基本统计分析方法6.1统计描述与区间估计6.6广义线性模型6.2均数比较的假设检验6.7生存分析6.3率比较的假设检验6.8混合效应模型6.4非参数假设检验6.9

练习案例6.5相关与回归6.10小结6.1统计描述与区间估计6.1统计描述与区间估计在做统计推断之前，一般需要对数据中的重要变量做统计描述。数值资料可以采用平均数、标准差等指标描述其集中位置或离散趋势，分类资料一般采用频数表来进行描述。在R语言中，置信区间（confidenceinterval）的估计一般由假设检验函数伴随给出（演示）。资料类型函数功能数值资料summary(x)极值、四分位数和平均数数值资料mean(x)算数均数数值资料median(x)中位数数值资料quantile(x,probs)百分位数数值资料sd(x)标准差数值资料var(x)方差数值资料IQR(x)四分位数间距数值资料min(x)最小值数值资料max(x)最大值分类资料table(x)根据N个变量创建N维频数表分类资料xtabs(formula,data)根据公式和数据创建N维频数表分类资料ftable(table)创建平铺式列联表表6.1常用统计描述函数6.2

均数比较的假设检验6.2均数比较的假设检验6.2.1t检验在研究中经常会对两组数据进行均数的比较（单因素，两水平），如临床试验中新药组的患者疗效是否优于传统药物组的药物疗效？患者服用降血压药物前的血压和服药后的血压是否发生了明显的变化？Student’st检验可以解决此类问题。在满足独立性、正态性、方差齐性的条件下，在R语言中可以使用t.test()进行t检验（演示）。如果想要比较当前样本与总体平均数间的差异，则可以使用单样本t检验（演示）。如果想要比较两组之间的差异，则可以使用成组t检验（two-samplet-test）进行平均数比较（演示）。若数据组内同质性差，则采用配对设计更为合适，只需在函数中添加paired=TRUE即可。6.2均数比较的假设检验6.2.1t检验类型格式说明单样本t检验t.test(x,mu)数据x与总体均数mu进行比较配对t检验t.test(x,y,paired=TRUE)数据x与数据y进行配对t检验成组t检验t.test(x,y,paired=FALSE)数据x与数据y进行成组t检验成组t检验t.test(y~x,data)数据集data中y在不同组别x中的成组t检验Shapiro-Wilks正态性检验shapiro.test(x)检验数据是否正态分布方差齐性检验bartlett.test(x,g,data)数据集data中y在不同组别x中的方差齐性表6.2R语言中的t检验实现方法6.2均数比较的假设检验6.2.2方差分析如果研究中的数据组数超过了2组，使用多次t检验则会造成一类错误（typeIerror）膨胀。方差分析（ANOVA）可以用于两个或两个以上水平的均数比较，甚至可以分析多个处理因素的效应，因此得到了广泛应用。在R语言中，方差分析可以通过aov()实现，基本语法为aov(formula,data)。符号用法~分隔符号，左边为观测值，右边为分组变量。如分析y在A、B两因素间的差异，写为y~A+B:表示交互项。如分析y在A、B两因素间及其交互作用的差异，写为y~A+B+A:B*表示所有可能交互项。y~A*B*C等价于y~A+B+C+A:B+B:C+A:C^交互项达到的次数。y~(A+B+C)^2等价于y~A+B+C+A:B+B:C+A:C.表示包含数据中除了y之外的所有变量。如果一个数据包含y,A,B,C,则y~.等价于y~A+B+C表6.4R语言中aov函数特殊符号6.2均数比较的假设检验6.2.2方差分析单因素方差分析：比较一个处理因素对应多个水平的总体均数差异（演示）。两因素方差分析：比较两个处理因素、多个水平的均数差异（演示）。若要分析两因素间的交互作用，可以使用*符号考虑交互项（演示）。表6.5R软件中不同类型方差分析实现方差分析类型表达式单因素方差分析aov(y~A)两因素方差分析aov(y~A+B)两因素方差分析（带交互作用）aov(y~A+B+A:B)重复测量方差分析,W为组内因子，B为组间因子aov(y~B*W+Error(Subject/W))6.2均数比较的假设检验6.2.2方差分析尽管方差分析可提示多组总体均数有显著差异，但并不能揭示究竟是哪些组间存在显著差异。此时可进行事后的两两比较，如Tukey'sHonestSignificantDifference（HSD）法（演示）。常用的两两比较的统计学方法还有SNK法和Dunnett法。SNK法可借助agricolae程序包中的SNK.test()实现（演示）。Dunnett法可借助multcomp程序包中的glht()实现（演示）。方差分析的另一种常见形式为重复测量方差分析，即受试者被测量了不止一次（演示）。6.3

率比较的假设检验6.3率比较的假设检验对于定性资料率的比较，卡方检验（Chi-SquareTest）是十分常用的方法，在R语言中可使用chisq.test()进行检验。趋势性检验可采用Cochran-Armitage法，借助DescTools程序包的CochranArmitageTest()；分层卡方检验可采用Cochran-Mantel-Haenszel法，在R语言中均可直接使用。若进行了多次卡方检验，想要进行事后两两比较，可直接采用假设检验多重校正法进行，如Benjamini&Hochberg法，使用p.adjust(P,method=“fdr”)；或者Bonferroni法，使用p.adjust(P,method=“bonferroni”)等。表6.6R语言中率比较的实现方法分析类型格式卡方检验chisq.test(table)Fisher确切概率法fisher.test(table)Cochran-Mantel-Haenszel检验mantelhaen.test(table)Cochran-Armitage趋势性检验CochranArmitageTest(table)6.3率比较的假设检验使用卡方检验需要满足一定的条件，对于四格表而言：（1）N＞40，且最小理论频数

＞5，直接使用卡方检验，即chisq.test(table,correct=F)（演示）；（2）N＞40，且最小理论频数≤5且

＞1，使用校正卡方检验，即chisq.test(table,correct=T)；（3）N≤40或最小理论频数≤1，使用Fisher确切概率法，即fisher.test(table)（演示）。进行Cochran-Mantel-Haenszel卡方检验可以使用mantelhaen.test()进行，用于检验两个名义变量在第三个变量的每一层中是否都是条件独立的。6.4

非参数假设检验6.4非参数假设检验经典的参数检验方法往往需要数据满足一定的前提假设，如t检验、方差分析需要数据满足独立性、正态性和方差齐性。如无法满足，非参数检验方法则是最合适的选择。非参数检验对于数据分布没有要求，适用范围极广，包括但不限于：（1）等级资料；（2）定量资料：正态分布但组间方差不齐、偏态分布、分布不明；（3）兼有定量与等级资料的性质；（4）资料中含有不确定值。但若该资料满足参数检验的条件，理应优先选择参数检验，否则会降低检验效能。6.4非参数假设检验若两组数据独立，属于完全随机设计，可以使用Wilcoxon秩和检验（或称为Mann-WhitneyU检验）来评估两组数据的总体分布是否相同（演示）。若两组数据不独立，属于配对设计，则可以使用Wilcoxon符号秩和检验来评估两组数据的总体中位数是否相同，只需在wilcox.test()中加入参数paired=TRUE即可（演示）。若数据组数多于两组，此时Wilcoxon秩和检验不再适用。如数据不满足方差分析的使用条件，则可以使用Kruskal-Wallis非参数检验来进行比较（演示）。Kruskal-Wallis检验类似于方差分析，显著差异的结论仅能表明多个组的分布不同或不全相同，若要进行事后两两比较，可使用Steel-Dwass检验（演示）。表6.6R语言中率比较的实现方法分析类型格式Wilcoxon秩和检验wilcox.test(y~x,data)Wilcoxon符号秩和检验wilcox.test(x,y,paired=TRUE)Kruskal-Wallis检验kruskal.test(y~x,data)6.5

相关与回归6.5相关与回归6.5.1相关相关（correlation）可以用于评估两变量间的共变关系。R语言可以计算多种相关系数，定量资料和等级资料可以采用cor()计算，方式为：cor(x,y,method)，对于定量资料一般采用Pearson相关分析，method指定为“pearson”即可，对于等级资料一般采用Spearman或Kendall秩相关分析，method指定为“spearman”或“kendall”即可。cor()仅可估计相关系数，优点是可以同时估计矩阵所有变量的两两相关系数；假设检验则可以使用cor.test(x,y,method)，但仅可检验两变量间的相关关系（演示）。若想要同时得到多个变量间的相关系数和假设检验结果，可借助psych程序包中的corr.test()。6.5相关与回归

6.5相关与回归6.5.2回归在R语言中，拟合线性回归模型的函数为lm()，格式为：fit<-lm(formula,data)（演示）。其中，formula是指要拟合模型的形式，data是一个数据框，包含用于拟合模型的数据。formula基本形式如下：

表6.9R语言中线性回归模型常用符号符号说明~分隔符号，~左边为因变量，~右边为自变量+分隔多个自变量:分析自变量间交互作用。如分析X1与X2间的交互，可写作Y~X1+X2+X1:X2*表示纳入所有可能的交互项6.5相关与回归6.5.2回归线性回归模型的方差分析表可以使用anova()查看。在很多情况下，需要分析的自变量数目往往不止一个，此时需要使用多元线性回归模型，依然使用lm()，即在~号的右边将想要分析的自变量使用+号连接起来。若想要分析全部潜在影响因素，则可使用fit<-lm(X~.,data)（演示）。表6.10线性回归模型具体信息展示函数说明summary()列出模型更具体的信息，如模型系数假设检验结果coefficients()列出模型的系数和截距项confint()列出模型系数和截距项的置信区间fitted()产生模型的预测值predict()根据新数据产生模型的预测值residuals()产生模型的残差anova()列出线性模型的方差分析表6.5相关与回归6.5.2回归若分析全部潜在影响因素的多元线性回归模型发现，没有任何一个因素具有统计学意义(p＞0.05)，很可能是由于部分自变量存在高度共线性，从而导致系数假设检验结果异常。此时，可进一步使用逐步回归（stepwiseregression）法进行变量筛选，使用step()实现（演示）。在线性回归模型中，任意变量之间一定存在着不同程度的相关性，这被称为多重共线性（Multicollinearity）。多重共线性是指线性回归模型中的解释变量之间存在精确相关关系或高度相关关系，从而使模型估计失真或难以估计准确。它会极大的影响线性回归模型估计结果的稳定性，导致系数估计有误差，甚至反常，模型结果解释起来非常困难。多重共线性可用统计量方差膨胀因子（VarianceInflationFactor，VIF）进行检测。VIF的平方根表示变量回归参数的置信区间能膨胀为与模型无关的预测变量的程度。car程序包中的vif()可计算VIF值。一般情况下，

VIF的平方根大于4表明存在多重共线性问题（演示）。6.6

广义线性模型6.6广义线性模型

6.6广义线性模型表6.11广义线性模型glm()相关参数模型类型分布族连接函数Logistic回归binomial(link="logit")多重线性回归gaussian(link="identity")Gamma回归Gamma(link="inverse")逆高斯回归inverse.gaussian(link="1/mu^2")泊松回归poisson(link="log")类高斯分布quasi(link="identity",variance="constant")类二项分布quasibinomial(link="logit")类泊松回归quasipoisson(link="log")6.6广义线性模型6.6.1Logistic回归评估若干个自变量与二分类型因变量的关联性，Logistic回归是一个十分常用的模型，使用glm()中的family=binomial()建立Logistic回归模型（演示）。此外，可结合使用exp()和coef()获取Logistic回归模型比值比（oddsratio，OR），使用confint()获取OR的置信区间（演示）。如果想要使用Logistic回归模型预测每个病例的具体复发概率，可以使用predict()，即可得到复发概率prob，函数格式为：bladder$prob<-predict(model,type="response")（演示）。6.6.2泊松回归当通过一系列连续型和/或分类型预测变量来预测计数型结果变量时，泊松回归（Poissonregression）是一个非常有用的工具。使用glm()中的family=poisson()建立Logistic回归模型（演示）。6.7

生存分析6.7生存分析生存分析（survivalanalysis）模型可用于处理“时间-事件”数据（time-to-eventdata），对个体一段时间内的风险进行预测，是临床研究中的常用方法。生存分析的显著特征是既考虑结局，又考虑生存时间，可充分利用截尾数据提供的不完全信息进行统计描述和统计推断。在R语言中，生存分析主要采用survival程序包进行。6.7生存分析6.7.1绘制生存曲线在survival程序包中，生存资料可以采用Surv()编码（演示）。描述生存数据最常用的方法是Kaplan-Meier生存曲线，可以使用survfit()绘制，若无分组，可以使用“~1”表示（演示）。可直接使用plot()绘制Kaplan-Meier生存曲线和95%置信区间（演示）。若有多个分组，可将分组变量添加到“~”右边，并可使用summary()查看指定时间点的事件数、生存率等信息（演示）。生存曲线的比较通常采用log-rank检验，可使用survdiff()进行（演示）。6.7生存分析

6.8

混合效应模型6.8混合效应模型混合效应模型（Mixedeffectmodel），或称多水平模型（Multilevelmodel），是现代回归分析的一种重要方法。传统回归需要满足的假设之一就是观测值之间是相互独立的，但很多情况下独立性并不能够被满足。例如，有一种实验设计称为重复测量设计（repeatedmeasuresdesign），或者称为被试内设计（within-subjectsdesign），这种设计方法中，同一受试对象会被测量多次，因为很多的测量都是来自同一受试对象，所以我们有理由相信这些测量之间是不独立的。诸如此类的不独立情况就不能使用传统回归模型，而混合效应模型就可以将不同水平（测量水平和个体水平）的变异都系统地分离开。6.8混合效应模型在R语言中，进行混合效应模型分析需要借助lme4或lmerTest等程序包，后者可同时输出变量的假设检验结果，因此更为推荐。以lmerTest程序包为例，线性混合效应模型函数为lmer(formula,data)；广义线性混合效应模型函数为glmer(formula,data)（演示）。建模方法如下：DV~IV+(1+R.Slope|R.Factor)其中，IV为固定因子，用于分析固定效应；R.Slope为随机因子对应的随机斜率；R.Factor为随机因子，用于分析随机效应；1|R.Factor为随机截距，表示不同随机因子对应的因变量分布不同。可使用anova()对固定效应进行检验，使用ranova()对随机效应进行检验（演示）。可以使用coef(summary())提取固定效应因子每个水平的假设检验结果（演示）。对于广义混合效应模型，glmer()用法与lmer()类似（演示）。6.9

练习案例6.9练习案例下面利用一个模拟数据集来探讨常见危险因素与肺癌发病之间的关联。该数据集有5000名研究对象，其中1000名来自肺癌病例组，4000名来自健康对照组。该数据集共包含12个变量，分别为研究对象的年龄、性别、BMI、民族（汉族/少数民族）、教育状况（高中及以上/高中以下/未知）、体育锻炼（经常/偶尔/从不）、是否吸烟、是否患COPD、是否患高血压、是否患糖尿病、是否患高脂血症，以及是否为肺癌新发病例等。数据结构如下，请以该数据集为例，结合本章R语言相关知识，选择合适的统计分析方法，探讨常见危险因素与肺癌发病之间的关联。编号年龄性别BMI高血压糖尿病高脂血症吸烟COPD民族教育状况体育锻炼肺癌162女性24是否是否是汉族高中及以上经常是228男性26否否是是否汉族高中及以上从不是···500080女性25否是否否否汉族高中及以上偶尔否表6.12数据结构6.10

小结6.10小结R语言是统计分析领域中备受欢迎的软件，其内置了丰富的经典统计方法函数，大量经典、优秀的外置程序包让其成为数据分析领域的翘楚。本章围绕经典统计分析框架，从统计描述和统计推断两方面出发，介绍了大量经典的统计分析方法，可供数据分析人员参考。6.10小结本章重要函数函数程序包作用主要用法summary()stats定量资料统计描述summary(object,...)table()base分类资料频数统计table(x)t.test()statst检验t.test(x,...)aov()stats方差分析aov(formula,data=NULL,projections=FALSE,qr=TRUE,contrasts=NULL,...)chisq.test()stats卡方检验chisq.test(x,y=NULL,correct=TRUE,p=rep(1/length(x),length(x)),rescale.p=FALSE,simulate.p.value=FALSE,B=2000)wilcox.test()stats秩和检验wilcox.test(x,...)cor()stats相关分析cor(x,y=NULL,use="everything",method=c("pearson","kendall","spearman"))cor.test()stats相关系数假设检验cor.test(x,...)lm()stats一般线性模型lm(formula,data,subset,weights,na.action,method="qr",mod