#计量基础-统计学R语言初步使用手册本手册为本人统计学复习时随意整理旨在为不熟悉统计学和R语言的同学做简单科普。它默认你不熟悉 R不熟悉统计学术语看到公式会紧张希望先“看懂”、再“会用”、最后“会解释”所以这份手册的写法不是“默认你已经懂”而是尽量按下面这个顺序展开先说这是什么再说为什么要学它再说它和前面学过的内容有什么关系最后才给公式和代码[!tip]这门课最重要的不是背一堆定义而是养成一个固定问法这个方法想解决什么问题它的输入是什么它的输出是什么图和数字该怎么解释0. 先建立全局地图0.1 这门课到底在学什么这门课不是在学“很多互不相关的统计方法”。它其实一直在围绕一个核心问题打转当我们同时观察很多变量时怎样理解它们之间的联合结构例如你现在不再只看一个人的身高而是同时看身高体重腰围血压肌肉量脂肪量这时问题就不再是“某一个变量是多少”而是这些变量之间是不是有共同变化规律能不能把很多变量压缩成少数几个“方向”如果有分组哪些方向最能区分组如果分成两大类变量这两类变量之间有没有整体关系这就是多元统计multivariate statistics的核心。0.2 这门课的方法之间是什么关系多元统计先理解数据怎么共同变化PCA压缩变量FA解释潜在结构LDA / QDA区分已知组别CCA研究两组变量的整体关系PLS / PLS-DA高维时更实用的替代FDR / q-value做很多检验时控制假阳性0.3 你真正需要达到的水平对期中来说你不需要做到能手推所有证明能背出特别复杂的推导能写出很长的 R 程序你需要做到的是知道每种方法在干什么能读懂基础输出能解释简单图形能比较相近方法的区别能在题目里选对方法1. R 语言最少必要入门这一章只讲本课立即会用到的 R。1.1 什么是 RR 是一种专门用于数据分析和统计的编程语言。你可以把它理解成一个计算器一个会画图的工具一个能做统计分析的软件语言对这门课来说R 的主要作用是组织数据调用统计函数画图帮助理解1.2 什么是“对象”在 R 里你做出的东西通常都可以存进一个名字里这个名字就叫对象object。x-1:5y-c(2,4,6,8,10)z-xy z这里x是一个对象y是一个对象z也是一个对象-的意思可以读成“把右边的结果放进左边这个名字里”1.3 什么是向量vector向量就是一串同类型的数据。例如一组身高一组分数一组年龄x-c(3,5,7,9)mean(x)sd(x)length(x)解释c(...)表示把数字合并成一个向量mean(x)求均值sd(x)求标准差length(x)求长度也就是里面有几个数1.4 什么是矩阵matrix矩阵就是一个长方形数字表。本课里你应该把矩阵理解成行 样本observations列 变量variablesX-matrix(c(160,55,165,60,170,68,175,74),nrow4,byrowTRUE)colnames(X)-c(height,weight)X dim(X)这段代码的意思有 4 个人每个人测了 2 个变量第 1 列是身高第 2 列是体重dim(X)会告诉你矩阵大小比如4 x 2。1.5 什么是数据框data.frame数据框比矩阵更灵活。矩阵通常要求所有列是同一种类型最常见是数值。数据框允许一列是数值一列是文字一列是类别标签df-data.frame(heightc(160,165,170,175),weightc(55,60,68,74),groupc(A,A,B,B))df这里height、weight是数值变量group是类别变量1.6 怎么取数据这一点很重要因为你后面几乎一直在做这个动作。df$height# 取名为 height 的这一列df[1,]# 第 1 行df[,2]# 第 2 列X[2,1]# 第 2 行第 1 列记忆方式[行, 列]前面空着 所有行后面空着 所有列1.7 本课最常用的基础函数mean(df$height)var(df$height)sd(df$height)cov(df$height,df$weight)cor(df$height,df$weight)scale(df[,c(height,weight)])这几个函数一定要看懂函数它在算什么你应该怎么理解mean()均值数据的平均水平var()方差数据波动有多大sd()标准差方差的平方根和原变量单位一致cov()协方差两个变量是不是一起变化cor()相关系数标准化后的协方差范围固定在 -1 到 1scale()标准化把变量转成可比较的统一尺度1.8 什么是apply()当你想对矩阵每一行或者每一列做同样的操作时apply()很方便。X-matrix(1:12,nrow3,byrowTRUE)apply(X,1,mean)# 每一行求均值apply(X,2,mean)# 每一列求均值这里的1和2很值得记住1 rows2 columns1.9 第一幅图散点图x-c(1,2,3,4,5)y-c(2,4,5,4,6)plot(x,y,pch19,colsteelblue,mainSimple scatterplot,xlabx,ylaby)这张图里每个点代表一个观测对象你要养成的读图习惯点大致沿一条斜线排列吗点是聚成几团吗有没有特别远的点后面 PCA、LDA、CCA 的很多图本质上都是“散点图升级版”。1.10 图示为什么标准化前后会长得不一样df-data.frame(incomec(40,50,55,60,70),heightc(155,160,165,170,175))par(mfrowc(1,2))plot(df$income,df$height,pch19,coltomato,mainRaw scale,xlabincome,ylabheight)z-scale(df)plot(z[,1],z[,2],pch19,colsteelblue,mainStandardized scale,xlabincome (z),ylabheight (z))par(mfrowc(1,1))左图和右图的区别不是“关系变了”而是“坐标尺度统一了”。这句话很重要标准化主要是为了让不同单位的变量能放在同一个公平的比较框架里。1.11 这一章结束后你应该会说什么什么是对象、向量、矩阵、数据框什么叫“行是样本列是变量”cov()和cor()都是在研究两个变量的关系scale()是把变量变到统一尺度上2. 最少必要统计基础这一章是后面所有方法的地基。如果这里没看懂后面的 PCA、CCA、FA 会显得像很多新名词如果这里看懂了你会发现后面大部分方法只是“换一种方式研究变量如何一起变化”。2.1 先分清几个最常见的统计对象变量variable变量就是你测量的某一项。例如身高体重年龄肌肉量观测observation观测就是一个样本、一位受试者、一条记录。例如第 1 个学生第 1 个病人第 1 个城市连续变量和分类变量连续变量可以取很多数值例如身高、体重、收入分类变量表示类别例如性别、物种、组别为什么这很重要因为方法选择会受变量类型影响连续变量常配 PCA、FA、CCA、PLS已知分类标签时常配 LDA、QDA、PLS-DA2.2 均值、方差、标准差均值是什么均值mean就是平均数。如果一组数是2, 4, 6那它们的均值就是(2 4 6) / 3 4均值回答的问题是这一组数据大概在什么位置方差是什么方差variance回答的问题是这些数离平均值有多分散如果所有数都非常接近均值方差就小。如果数值很分散方差就大。样本方差公式s 2 ∑ ( x i − x ˉ ) 2 n − 1 s^2 \frac{\sum (x_i - \bar{x})^2}{n - 1}s2n−1∑(xi​−xˉ)2​这公式你不需要一上来死背但你要知道每一部分在干什么x_i - xbar每个数离平均值有多远平方避免正负抵消求和把所有偏离加起来除以n - 1得到样本方差标准差是什么标准差standard deviation就是方差的平方根。它比方差更常被直观理解因为标准差和原数据单位一致例如身高的标准差还是“厘米”收入的标准差还是“元”最小例子x-c(2,4,6)mean(x)var(x)sd(x)2.3 协方差两个变量是否一起变化先别看公式先看意思协方差covariance想回答的是当变量 X 变大时变量 Y 会不会也跟着变大如果是就倾向于正协方差。如果 X 变大时 Y 倾向于变小就倾向于负协方差。公式长什么样cov ⁡ ( X , Y ) ∑ ( x i − x ˉ ) ( y i − y ˉ ) n − 1 \operatorname{cov}(X, Y) \frac{\sum (x_i - \bar{x})(y_i - \bar{y})}{n - 1}cov(X,Y)n−1∑(xi​−xˉ)(yi​−yˉ​)​不要被公式吓到它本质上仍然是在看X 离自己均值的偏差Y 离自己均值的偏差这两个偏差是不是经常同号如果常常同号一个高于平均时另一个也高于平均协方差倾向于正如果常常异号一个高于平均时另一个低于平均协方差倾向于负例子x-c(1,2,3,4)y-c(2,4,6,8)z-c(8,6,4,2)cov(x,y)cov(x,z)解释cov(x, y)应该是正的因为它们一起增大cov(x, z)应该是负的因为一个变大时另一个变小2.4 相关更容易解释的关系强度为什么还需要相关因为协方差有一个问题它会受单位影响例如你的收入用“元”记录另一个人用“万元”记录同样的关系协方差数值可能完全不同。所以我们引入相关系数correlation。相关的好处相关系数的范围固定在-1到1这样解释就容易得多1完全正线性关系-1完全负线性关系0没有明显线性关系最小例子cor(x,y)cor(x,z)2.5 图示高相关和低相关看起来有什么不同set.seed(767)x-rnorm(100)y1-0.9*xrnorm(100,sd0.2)y2-rnorm(100)par(mfrowc(1,2))plot(x,y1,pch19,colsteelblue,mainStrong linear relation)plot(x,y2,pch19,coltomato,mainWeak linear relation)par(mfrowc(1,1))cor(x,y1)cor(x,y2)这段代码最应该帮你建立的是“图形直觉”左图里点大致沿斜线分布相关高右图里点更像随机一团相关低这件事很重要因为后面很多方法其实都在研究这种“联合变化的形状”。2.6 标准化为什么重要这是很多初学者第一次真正卡住的地方所以这里讲慢一点。第一步先理解“标准化”这个词标准化standardization的意思不是“让数据变得标准答案一样”。它真正做的事情是把不同单位、不同量纲的变量转换到同一个可比较的尺度上。标准化公式是z x − x ˉ s d ( x ) z \frac{x - \bar{x}}{sd(x)}zsd(x)x−xˉ​你可以把这条公式拆成两步理解x - xbar先看这个值离平均值有多远再除以sd(x)把“远多少”转成“相对于本变量自己的波动来说远多少”第二步标准化后会发生什么标准化之后均值会变成 0标准差会变成 1这不代表数据“失真”了而是代表它们被放进了同一种单位体系里第三步为什么 PCA 特别在意这件事先想一个最直观的例子。你现在有两个变量income收入数值大概是40, 50, 60height身高数值大概是155, 165, 175如果你直接把它们放进 PCA不标准化会怎样答案是数值范围更大的变量更容易主导分析这时候 PCA 很可能更关注“数字大不大”而不是“变量结构是否重要”。也就是说不标准化时结果可能更多反映变量单位而不是变量之间真正的关系第四步最小代码例子df-data.frame(incomec(40,50,55,60,70),heightc(155,160,165,170,175))scale(df)你运行后会看到原来的income和height数值都变了但它们都被改造成“均值 0、标准差 1”的尺度第五步图形理解df-data.frame(incomec(40,50,55,60,70),heightc(155,160,165,170,175))par(mfrowc(1,2))plot(df$income,df$height,pch19,coltomato,mainBefore standardization,xlabincome,ylabheight)z-scale(df)plot(z[,1],z[,2],pch19,colsteelblue,mainAfter standardization,xlabincome (z),ylabheight (z))par(mfrowc(1,1))这段图的重点不是让你看“点云完全变形了没有”而是让你明白标准化以后两列变量终于可以在同一个公平框架里比较第六步什么时候特别要考虑标准化这些场景通常应该优先考虑变量单位不同变量数值范围差很多你不希望“大数字变量”天然拥有更大话语权第七步你现在应该会怎么说标准化就是把变量转换到同一个可比较的尺度上。它特别重要因为像 PCA 这类方法会受变量尺度影响如果不同变量的单位差很多不标准化时大量纲变量可能主导结果。2.7p-value的最小正确理解这一小节也常卡人因为很多教材上来就讲得很抽象。我们先只抓最小正确含义。第一步先理解零假设H0零假设H0可以先粗略理解为“没有差异”“没有关系”“没有效应”比如两组均值没有差异两个变量没有关系这不一定是现实真相而是检验开始时先假设它成立。第二步p-value到底在说什么p-value可以理解为如果零假设真的成立那么像现在这么极端的结果出现的概率有多大。注意它说的是“假设 H0 成立时数据有多极端”它没有说H0为真的概率是多少这两句话完全不是一回事。第三步一个硬币的直觉例子假设你怀疑一枚硬币不公平。零假设这枚硬币是公平的你抛了 20 次结果出现了 19 次正面。这时你会想如果硬币真是公平的出现这么偏的结果是不是太夸张了这时对应的p-value就是在回答在“硬币公平”的前提下像这么极端的结果有多罕见第四步p-value小意味着什么如果p-value很小通常说明在H0为真的前提下这个结果不太常见因此我们会倾向于怀疑H0但不要说成“H0 有 3% 的概率是真的”这是错误说法。第五步不要把p-value当成什么不要把它理解为零假设为真的概率效应大小结果是否重要的全部标准例如很小的效应如果样本非常大也可能有很小的p-value很大的效应如果样本很小也不一定有很小的p-value第六步你现在应该会怎么说p-value表示如果零假设为真那么观察到当前这么极端或更极端结果的概率有多大。它不是零假设为真的概率也不是效应大小。2.8 为什么会有多重检验问题这件事如果不靠例子讲会很容易抽象。第一步先看单个检验假设你做 1 次检验并规定p 0.05就算显著这大致意味着即使零假设为真也有 5% 左右概率因为随机波动而误报“显著”这个错误就叫假阳性false positive。第二步如果你做很多很多次检验现在不是做 1 次而是做 1000 次。就算这 1000 个检验全都没有真实效应也会发生什么答案是你仍然可能看到很多“显著”为什么因为每次检验都有机会因为随机波动误报。第三步最直观的数字感觉如果每次检验的误报概率大约是 5%那做 1000 次检验时你大概可能看到几十个看起来“显著”的结果虽然它们其实只是随机误报。这就是为什么检验做得越多只看原始p-value就越危险第四步为什么这在高维数据里特别严重在多元统计和现代数据分析里经常会出现几百个变量几千个基因很多后续单变量检验这时如果你还按“每个检验都单独看0.05”的方式判断就很容易把噪声当成信号。第五步这就是 FDR 为什么重要FDRFalse Discovery Rate想控制的是你判为显著的一堆结果里有多少比例可能是假阳性它不是在说“完全不许犯错”而是在说我们要控制错误发现的比例不要让假阳性泛滥第六步最小模拟例子set.seed(767)pvals-runif(1000)sum(pvals0.05)这里runif(1000)生成了 1000 个 0 到 1 之间的随机数。你可以把它粗略看成在零假设都成立时很多检验产生的一组原始p-valuesum(pvals 0.05)统计的是其中有多少个“看起来显著”你通常会发现明明没有真实信号仍然会冒出一些“显著”这就是多重检验问题最直观的样子。第七步你现在应该会怎么说当同时进行很多检验时即使每个检验都没有真实效应也可能因为随机波动出现不少看似显著的结果所以不能只看单个原始p-value这就是 FDR 等多重检验控制方法存在的原因。2.9 这一章结束后你应该会说什么方差描述一个变量的波动大小协方差描述两个变量是否一起变化相关是更容易解释的标准化关系强度标准化是把变量放到同一个可比较尺度上p-value不是零假设为真的概率检验做很多次时会自然冒出假阳性所以要控制多重检验错误