Size Matters: Metabolic Rate and Longevity (Regression analysis sample)

Size Matters: Metabolic Rate and Longevity

John Tukey once said, “The best thing about being a statistician is that you get to play in everyone’s backyard.” I enthusiastically agree!

I frequently enjoy reading and watching science-related material. This invariably raises questions, involving other “backyards,” that I can better understand using statistics. For instance, see my […]

数据分析之美:如何进行回归分析

1. 确定自变量与Y是否相关 证明:自变量X1,X2,….XP中至少存在一个自变量与因变量Y相关 For any given value of n(观测数据的数目) and p(自变量X的数目), any statistical software package can be used to compute the p-value associated with the F-statistic using this distribution. Based on this p-value, we can determine whether or not to reject H0. (用软件计算出的与F-statistic 相关的p-value来验证假设,the p-value associated with the F-statistic) 例子: Is there a relationship between […]

STAMP (Statistical Analysis of Metagenomic Profiles)

STAMP

Using STAMP to identify SEED subsystems which are differentially abundant between CandidatusAccumulibacter phosphatis sequences obtained from a pair of enhanced biological phosphorus removal (EBPR) sludge metagenomes(data originally described in Parks and Beiko, 2010).

STAMP (Statistical Analysis of Metagenomic Profiles) is a software package for analyzing metagenomic profiles (e.g., a taxonomic profile indicating the […]

Meta分析(Meta-analysis)简介

Meta分析的由来: Meta分析(Meta-analysis)是心理学家Glass于1976年首次提出,原文是这样的: Meta-analysis refers to the analysis of analyses. I use it to refer to the statistical analysis of a large collection of results from individual studies for the purpose of integrating the findings. It connotes a rigorous alternative to the casual, narrative discussions of research studies which typify our attempts to make sense of […]

真理在缩水,还是上帝在掷骰子?(转贴)

最近在Google Reader中看见科学松鼠会有两篇文章被频繁分享,名为《真理在缩水——现代科学研究方法并不尽善尽美?》(上)与(下),下文简称《缩水》。文章很有意思,而实际上说的是我们的老本行——统计学,因此我在这里也发表一些我的想法和理解,包括这两年我在美帝学习的一些思考,部分内容受益于两位老师Kaiser和Nettleton教授,先向他们致谢(尽管他们永远都不会看到这篇文章)。同时我也要先说明一下,读这篇文章可能会很花时间(至少我花了大约二十小时写这篇文章),即使我的观点没有价值,我相信里面的引用文献是有价值的。

初读文章,我脑子里冒出的一句话是“上帝在跟我们掷骰子”,文中给出了大量的不可重复的试验,仿佛就像那些号称“具有统计学意义”(下文我再说这个所谓的“意义”)的试验结果只是若干次骰子中的一次运气好的结果而已。读完文章,我们可能不禁要问,到底是真理在缩水,还是它根本就不曾存在?下面我从四个方面来展开,分别说明人对随机性的认识、统计推断的基石、让无数英雄折腰的P值、以及可重复的统计研究。

一、感知随机

随机变量在统计分析中占据中心地位,数学上关于随机变量的定义只是一个“干巴巴的函数”,从样本空间映射到实数集,保证从实数集上的Borel域逆回去的集合仍然在原来的sigma域中即可。随机变量的性质由其分布函数刻画。写这段话的目的不是为了吓唬你,也不是为了作八股文,而是来说明我为什么不喜欢数学的理由,对我而言,我觉得有些数学工具只是为了让自己不要太心虚,相信某时某刻某个角落有个理论在支撑你,但后果是弱化了人的感知,当然,也有很多数学工具有很强的直觉性(如果可能,我想在未来下一篇文章里面总结这些问题)。我一直认为很多人对随机性的感知是有偏差的,对概率的解释也容易掉进陷阱(参见Casella & Berger的Statistical Inference第一章,例如条件概率的三囚徒问题)。

《缩水》一文发表了很多不可重复的试验案例,我们应该吃惊吗?我的回答是,未必。举两个简单的例子:

第一个例子:很多数据分析人员都很在意所谓的“离群点”,论坛上也隔三差五有人问到这样的问题(如何判断和处理离群点),而且也有很多人的做法就是粗暴地删掉,我从来都反对这种做法。除了基于“数据是宝贵的”这样简单的想法之外,另一点原因是,离群点也许并非“异类”。离群点是否真的不容易出现?请打开R或其它统计软件,生成30个标准正态分布N(0, 1)随机数看看结果,比如R中输入rnorm(30),这是我运行一次的结果:

> rnorm(30) [1] 1.19062761 -0.85917341 2.90110515 0.59532402 -0.05081508 -0.06814796 [7] 2.08899701 0.76423007 0.92587075 -1.16232929 -0.68074378 -1.40437532 [13] -0.17932604 -0.72980545 -0.53850923 0.21685537 -0.35650714 -1.32591808 [19] -0.88071526 -1.25832441 0.24001498 -0.41682799 -0.09576492 -0.17059052 [25] -0.99947485 0.25108253 -0.47566842 -0.28028786 0.79856649 -0.13250974

30在现实中是一个比较小的样本量,你看到了什么?2.901?它接近3倍标准差的位置了。还有2.089?……如果你不知道这批数据真的是从标准正态分布中生成出来的,现在你会有什么反应?把2.9删掉?标准正态分布是一个在我们眼中很“正常”的分布,而一个不太大的样本量一次试验足以生成几个“离群点”,那么要是成千上万的试验中没能产生几项震惊世界的结果,你会怎样想?(上帝的骰子坏掉了)

另一个例子和统计学结合紧密一点,我们谈谈残差的QQ图。QQ图是用来检查数据正态性的一种统计图形,与腾讯无关,细节此处略去,大意是图中的点若呈直线状(大致分布在对角线上),那么可以说明数据的正态性比较好,因此QQ图经常被用在对回归模型残差的正态性诊断上。我的问题是,即使数据真的是正态分布,你是否真的会看见一些分布在直线上的点?若答案是否定的,那么我们就得重新审视我们对分布和随机的认识了。下图是一幅教科书式的QQ图(仍然基于30个正态分布随机数):

“正常的”QQ图(来自R语言qqnorm(rnorm(30)))

随机性并没有这么美好,即使数据真的来自正态分布,你也有可能很容易观察到歪歪扭扭的QQ图,尤其是小样本的情况下。比如下图是50次重复抽样的正态数据QQ图,它和你想象的QQ图本来的样子差多远?

library(animation) set.seed(710) […]

Using prcomp/princomp for PCA in R (三)

Testing i.pca ~ prcomp(), m.eigensystem ~ princomp()

1. Briefly about PCA 2. The modules/functions that implement PCA in GRASS & R 3. My claims (Entitled Comments) 4. Evidence (=the numbers derived from i.pca, prcomp, princomp, m.eigensystem using some MODIS surface reflectance bands).

Finally all is clear _but_ one thing: the only “unknown” variable (to […]

Using prcomp/princomp for PCA in R (二)

###############################

PCA ############################### install.packages(“vegan”) library(vegan)

> STpcoa<-read.table(file=”bactera_16s_final.subsample.phylip.tre1.weighted.phylip.pcoa.axes”, header=T,row.names=1) > STpcoa axis1 axis2 axis3 axis4 Cellulose -0.020878 -0.234601 0.167454 0 Foodwaste -0.234592 0.221741 0.085802 0 Sludge 0.368882 0.100725 -0.010570 0 Xylan -0.113413 -0.087865 -0.242686 0 >pl.STpcoa<-princomp(STpcoa) > summary(pl.STpcoa) Importance of components: Comp.1 Comp.2 Comp.3 Comp.4 Standard deviation 0.2260563 0.1746944 0.1536319 0 Proportion of Variance 0.4856521 0.2900347 […]

Using prcomp/princomp for PCA in R (一)

Difference between prcomp and princomp:

‘princomp’ can only be used with more units than variables”

prcomp是基于SVD分解(svd()函数,princomp是基于特征向量eigen()函数)

Good video source:

http://www.youtube.com/watch?v=oZ2nfIPdvjY

http://www.youtube.com/watch?v=I5GxNzKLIoU&feature=relmfu

http://www.planta.cn/forum/viewtopic.php?t=16754&highlight=%D3%EF%D1%D4

###########################################

以下所有代码包括练习数据,都可在R平台上直接运行。

#主成分分析和主成分回归 主成分分析的思想是Pearson 1901年提出的,Hotelling 1933进一步发展 在R中,进行主成分分析用到princomp() 函数

用法 princomp(x, cor = FALSE, scores = TRUE, covmat = NULL, subset = rep(TRUE, nrow(as.matrix(x))), …)

# 分析用数据 # cor 是否用样本的协方差矩阵作主成分分析 prcomp() 二 summary()函数 三 […]

4sample CA RDA analysis

 

> gtsdata_test=read.table(“gtsdata.txt”, header=T) > gtsenv=read.table(“gtsenv.txt”, header=T) > gtsdata_data_t<-t(gtsdata_data) > decorana(gtsdata_data_t)

Call: decorana(veg = gtsdata_data_t)

Detrended correspondence analysis with 26 segments. Rescaling of axes with 4 iterations.

DCA1 DCA2 DCA3 DCA4 Eigenvalues 0.8634 0.4834 0.23788 0 Decorana values 0.8721 0.3793 0.07223 0 Axis lengths 5.3292 2.1115 1.80907 0

> gts.ca=cca(gtsdata_data_t) > gts.ca Call: cca(X = […]

基于Vegan 软件包的生态学数据排序分析学习

“基于Vegan 软件包的生态学数据排序分析 赖江山 米湘成 (中国科学院植物研究所植被与环境变化国家重点实验室,北京 100093) 摘要:群落学数据一般是多维数据,例如物种属性或环境因子的属性。多元统计分析是群落生态学常用的分析方法,排序(ordination)是多元统计最常用的方法之一。CANOCO是广泛使用的排序软件,但缺点是商业软件价格不菲,版本更新速度也很慢。近年来,R语言以其灵活、开放、易于掌握、免费等诸多优点,在生态学和生物多样性研究领域迅速赢得广大研究人员的青睐。R语言中的外在软件包“Vegan”是专门用于群落生态学分析的工具。Vegan能够提供所有基本的排序方法,同时具有生成精美排序图的功能,版本更新很快。我们认为Vegan包完全可以取代CANOCO,成为今后排序分析的首选统计工具。本文首先简述排序的原理和类型,然后介绍Vegan的基本信息和下载安装过程,最后以古田山24公顷样地内随机抽取40个20m×20m的样方为例,展示Vegan包内各种常用排序方法(PCA,RDA,CA和CCA)和排序图生成过程,希望能为R的初学者尽快熟悉并利用Vegan包进行排序分析提供参考。

gtsdata

gtsenv.txt

赖江山.pdf

> setwd(“/winxp_disk2/shenzy/R/Vegan”) > gtsdata=read.table(“gtsdata.txt”, header=T) > gtsenv=read.table(“gtsenv.txt”, header=T) > install.packages(“vegan”) Installing package(s) into ‘/home/shenzy/R/x86_64-pc-linux-gnu-library/2.15’ (as ‘lib’ is unspecified) 试开URL’http://cran.csiro.au/src/contrib/vegan_2.0-4.tar.gz’ Content type ‘application/x-gzip’ length 1576584 bytes (1.5 Mb) 打开了URL ================================================== downloaded 1.5 Mb * installing *source* package ‘vegan’ … ** 成功将‘vegan’程序包解包并MD5和检查 ** libs gfortran -fpic -O3 […]