summary()
sapply(x,fun,options):对数据框或矩阵中的每一个向量进行统计
mean
sd：标准差
var：方差
min:
max:
median:
length:
range:
quantile:
vars <- c("mpg", "hp", "wt")
head(mtcars[vars])
summary(mtcars[vars])
mystats <- function(x, na.omit = FALSE) {
if (na.omit)
x <- x[!is.na(x)] #去掉缺失值
m <- mean(x)
n <- length(x)
s <- sd(x)
skew <- sum((x - m)^3/s^3)/n #偏度
kurt <- sum((x - m)^4/s^4)/n - 3 #峰度
return(c(n = n, mean = m, stdev = s, skew = skew, kurtosis = kurt))
}
sapply(mtcars[vars], mystats)
#如结果，分布右偏(0.61)，较正态分布稍平(-0.37)
#Hmisc包中的describe可返回变量和观测的数量，缺失值和唯一值的数目、平均值、分位数及五个最大的值和最小的值
library(Hmisc)
describe(mtcars[vars])
#pastects包中有一个stat.desc函数，可以计算各种描述性统计量
library(pastecs)
stat.desc(mtcars[vars])
#psych包的describe可以计算描述性统计量
library(psych)
describe(mtcars[vars])
#aggregate()：分组获取描述性统计量
aggregate(mtcars[vars], by = list(am = mtcars$am), mean)
aggregate(mtcars[vars], by = list(am = mtcars$am), sd)
#by(data,indices,fun)
dstats <- function(x)(c(mean=mean(x), sd=sd(x)))
by(mtcars[vars], mtcars$am, dstats)
#doBy的summaryBy可以分组计算描述性统计量
library(doBy)
summaryBy(mpg + hp + wt ~ am, data = mtcars, FUN = mystats)
#psych中的desribe.by可以分组计算描述性统计量，但不允许指定任何函数。
library(psych)
describe.by(mtcars[vars], mtcars$am)
#使用reshape包
library(reshape)
dstats <- function(x) (c(n = length(x), mean = mean(x), sd = sd(x)))
#id.vars，按什么分组，measure看什么变量？
dfm <- melt(mtcars, measure.vars = c("mpg", "hp", "wt"), id.vars = c("am", "cyl"))
cast(dfm, am + cyl + variable ~ ., dstats)
table：使用N个类别型变量创建一个N维列联表
xtabs：根据一个公式和一个矩阵或数据框创建一个N维列联表
prop.table：将表中条目表示为分数形式
margin.table：计算表中条目的和
addmargins：将结果放入表中
ftable：创建一个紧凑的平铺式列联表
#一维列联表
mytable <- with(Arthritis, table(Improved))
mytable
prop.table(mytable)
#二维列联表
mytable <- xtabs(~ Treatment+Improved, data=Arthritis)
mytable
#计算边缘频数，1表示第一个变量
margin.table(mytable, 1)
#计算边缘比例，1表示第一个变量
prop.table(mytable, 1)
margin.table(mytable, 2)
prop.table(mytable, 2)
prop.table(mytable)
#每行和每列最后加sum
addmargins(mytable)
admargins(prop.table(mytable))
#仅添加各行的和
addmargins(prop.table(mytable, 1), 2)
#仅添加各列的和
addmargins(prop.table(mytable, 2, 1)
#可以用gmodels的CrossTable生成二维列联表
library(gmodels)
CrossTable(Arthritis$Treatment, Arthritis$Improved)
#独立性检验
library(vcd)
mytable <- xtabs(~Treatment+Improved, data=Arthritis)
chisq.test(mytable)
#以上代码结果p值为0.001，存在某种关系
mytable <- xtabs(~Improved+Sex, data=Arthritis)
chisq.test(mytable)
#以上代码结果p值为0.08，没关系
#以下代码结果p值为0.001，有关系
mytable <- xtabs(~Treatment+Improved, data=Arthritis)
fisher.test(mytable)
#以下代码p值为0.0006，假设不存在三阶交互作用
mytable <- xtabs(~Treatment+Improved+Sex, data=Arthritis)
mantelhaen.test(mytable)
#以下为相关性度量，值越大相关性越强
library(vcd)
mytable <- xtabs(~Treatment+Improved, data=Arthritis)
assocstats(mytable)
#cor函数：计算Pearson相关系数，Spearman相关系数，Kendall相关系数
#cov函数：计算协方差
states <- state.x77[, 1:6]
#以下代码计算方差和协方差
cov(states)
#以下代码计算Pearson积差相关系数
cor(states)
#以下代码计算Spearman等级相关系数
cor(states, method="spearman")
#可以看出收入和高中毕业率之间存在强正相关，文盲率和预期寿命之间存在很强的负相关
#设定行列，可以得出更直观的结果
x <- states[, c("Population", "Income", "Illiteracy", "HS Grad")]
y <- states[, c("Life Exp", "Murder")]
cor(x, y)
#相关性检验
cor.test(states[, 3], states[, 5])
#一次检验多个
library(psych)
corr.test(states, use = "complete")
#独立样本t检验
library(MASS)
t.test(Prob ~ So, data=UScrime).
#非独立样本t检验
library(MASS)
sapply(UScrime[c("U1", "U2")], function(x) (c(mean = mean(x),
sd = sd(x))))
with(UScrime, t.test(U1, U2, paired = TRUE))
#组间差异的非参数检验，wilcox.test适用于两组成对数据和无法保证正态性假设的情境
with(UScrime, by(Prob, So, median))
wilcox.test(Prob ~ So, data=UScrime)
#paired设为TRUE，配对检验
sapply(UScrime[c("U1", "U2")], median)
with(UScrime, wilcox.test(U1, U2, paired = TRUE))
#多于两组的检验
states <- as.data.frame(cbind(state.region, state.x77))
kruskal.test(Illiteracy ~ state.region, data=states)