R语言_基本统计分析

#基本统计分析

#整体描述性统计分析,针对数值变量
attach(mtcars)
opar = par(no.readnoly=TRUE)
d = mtcars[c("mpg","hp","wt")]
head(d)
#summary
#较标准正态分布呈现正偏,且较平。(偏度为正,峰度为负)
summary(d)
plot(density(mpg))
#describe
#多了峰度,偏度等数据
library(psych)
describe(d)


#分组描述统计,针对数值变量
#aggregate,fun = c(mean,sd)
aggregate(d,by=list(am=am),mean)
#by, func = self-designed
dstats <- function(x){
    c(mean=sapply(x,mean),sum=sapply(x,sum))
}
tapply(d$mpg,factor(am),dstats)
by(d,factor(am),dstats)
by(d,factor(am),summary)
#summaryBy
library(doBy)
dstats <- function(x){
    c(mean=mean(x),sum=sum(x))
}
summaryBy(mpg+hp+wt~am,data=mtcars,FUN=dstats)
#describe.by 
#多了峰度,偏度等数据,不允许指定函数,适应度较差
library(psych)
describeBy(d,am)
#melt and casr
library(reshape)
dstats <- function(x){
    c(mean=mean(x),sum=sum(x),length=length(x))
}
dfm = melt(mtcars,measure.vars = c("mpg","hp","wt"),
           id.vars=c("am","cyl"))
cast(dfm,am+cyl+variable~.,dstats)



#频数表和列联表,针对类别变量
#函数总概
table(var1,var2)
xtabs(formula,data) #根据一个公式和一个矩阵或者数据框创建n维列联表
prop.table(table,margins) #根据margins定义的边际列表将表中条目表示为分数形式
margin.table(table,margin) #依据margin定义的边界计算和
addmargins(table,margins) #将margin(默认求和结果)放入表中
ftable(table) #创建一个紧凑的平铺式的列联表
#一维列联表
#table默认忽略缺失值,若不则useNA="ifany"
table(cyl)
prop.table(table(cyl))
#二维列联表
table(cyl,am)
mytable = xtabs(~am+cyl,data=mtcars)
margin.table(mytable,margin = 1) #行
prop.table(mytable,margin = 1) #行
prop.table(mytable,margin = 2) #列
prop.table(mytable) #行列所占比例
#添加边际和的二维列联表
addmargins(mytable)
addmargins(prop.table(mytable))
addmargins(prop.table(mytable,1),2)
addmargins(prop.table(mytable,2),1)
#CrossTable生成二维列联表
install.packages("gmodels")
library(gmodels)
CrossTable(am,cyl)
#多维列联表
mytable = xtabs(~am+cyl+gear,data=mtcars)
mytable
ftable(mytable)
margin.table(mytable,1)
margin.table(mytable,3)
margin.table(mytable,c(1,3))
ftable(prop.table(mytable,c(1,2))) 
ftable(addmargins(prop.table(mytable,c(1,2)),3) )
#将表转化为扁平格式
table2flat <- function(mytable) {
    df = as.data.frame(mytable)
    rows = dim(df)[1]
    cols = dim(df)[2]
    x = NULL
    counts=0
    for (i in 1:rows) {
        for (j in 1:df$freq[i]) {
            counts = counts+1
            row = df[i,c(1:(cols-1))]
            x = rbind(x,row)
        }
        print(c(i,counts))
    }
    row.names(x) <- c(1:dim(x)[1])
    return (list(x,counts))
}

treatment = rep(c("Placebo","Treated"),times=3)
improved  = rep(c("None","Some","Marked"),each=2)
freq = c(29,13,7,17,7,21)
mytable = as.data.frame(cbind(treatment,improved,freq))
mytable$freq = as.numeric(as.character(mytable$freq))
mydata = table2flat(mytable)



#独立性检验,描述类别变量独立性
#卡方独立性检验
#卡方备注:
#p值表示从总体中抽取样本行变量与列变量相互独立的概率,
# p<0.01,概率非常小,所以拒绝相互独立的原假设
# p>0.05,概率不够小,没有足够理由说明原来的两个变量是不独立的
#产生警告的原因,是6个单元格(男性,一定程度改善)有一个小于5,可能使卡方无效
library(vcd)
mytable = xtabs(~Treatment+Improved,data=Arthritis)
chisq.test(mytable)  #治疗和改善不独立 p<0.01
mytable = xtabs(~Sex+Improved,data=Arthritis)
chisq.test(mytable)  #性别和改善独立 p>0.05
#Fisher精确检验
#原假设是:边界固定的列联表中行和列是相互独立的
mytable = xtabs(~Treatment+Improved,data=Arthritis)
fisher.test(mytable)
#Cochran-Mantel-Haenszel检验
#原假设是:两个名义变量在第三个变量的每一层中都是条件独立的
#下面检验治疗情况和改善情况在性别的每一个水平下是否独立,检验不存在三阶交互作用
#结果表明:患者接受的治疗与得到的改善在性别的每一个水平下并不独立
#即:分性别看,用药治疗的患者较接受安慰剂的患者得到更多改善
mytable = xtabs(~Treatment+Improved+Sex,data=Arthritis)
mantelhaen.test(mytable)
#相关性的度量
#上面的显著性评估目的是判断变量间是否相互独立
#若不,则接着衡量相关性的强弱
#共得到了phi,列联,Cramer‘s V系数,较大意味着相关性越强
library(vcd)
mytable = xtabs(~Treatment+Improved,data=Arthritis)
chisq.test(mytable)
assocstats(mytable)



#相关性
#上述的独立性检验主要描述类别变量的独立性
#针对定量变量,使用相关性去描述
#Pearson 衡量两个定量变量之间的线性相关程度
#Spearman 衡量分级定序变量之间的相关程度
#Kendall 非参数的等级相关度量
#cor(x,use=,method=) use制定缺失数据的处理方式,method制定相关系数的类型,默认use=everything,method=person
#但是结果并不指明相关系数是否显著不为0,即根据样本数据是否有足够的证据得出整体相关系数不为0
states=state.x77[,1:6]
cov(states)
cor(states)
cor(states,method="spearman")
x = states[,c("Population","Income","HS Grad")]
y = states[,c("Life Exp","Murder")]
cor(x,y)
#偏相关
#指控制一个或多个定量变量时,另外两个定量变量之间的相互关系
install.packages("ggm")
library(ggm)
pcor(c(1,5,2,3,6),cov(states))  #控制2,3,6的影响,判断1,5的相关系数
#相关性的显著检验
#原假设:变量不相关,相关系数为0
#cor.test(x,y,alternative=,method=)
cor.test(states[,3],states[,5])
#计算相关矩阵并进行显著性检验
library(psych)
corr.test(states,use="complete")



#t检验
#关注连续型变量的组间比较,类别型变量参考上文独立性检验部分
#例子:新药治疗的患者相比旧药是否有更大程度改善;新工艺是否比旧工艺制造的不合格产品更少
#独立样本的t检验
#假设:两个总体的均值相等,并且从正态总体中取得
#下面进行假设方差不等的双侧检验,比较南方和非南方的监禁概率
#可以拒绝相同监禁概率的假设,因为p<0.01
library(MASS)
t.test(Prob~So,data=UScrime)
#非独立样本的t检验
#假设:组件的差异呈现正态分布
#P值反映了如果实际均值相等,那么获得一个差异如此大的样本的概率小于2.2e-16
library(MASS)
sapply(UScrime[c("U1","U2")],function(x) (c(mean=mean(x),sd=sd(x))) )
with(UScrime, t.test(U1,U2,paired=TRUE))
#多于两组的情况
#假设数据从正态总体中独立抽样而得
ANOVA分析


#组件差异的非参数检验
#如果数据无法满足t检验或者anova的参数假设,一般采用非参数方法
#例如:结果变量在本质上就严重偏斜或呈现有序关系
#两组的比较
#若两组数据独立,可以使用Wolcoxon秩和检验(Mann-Whitney U检验)。来评估观测是否是从相同概率分布中抽的
#即:在一个总体中获得更高得分的概率是否比另一个总体更大
#评价:是非独立样本t检验的一种非参数替代方法。适用于两组成对数据和无法保证正态性假设的情景。
#1
with(UScrime,by(Prob,So,median))
wilcox.test(Prob~So,data=UScrime)
#2
#在本例中,含参的t检验和其作用相同的非参数检验得到了相同的结论。
#当t检验的假设合理时,参数检验的功效更强(更容易发现存在的差异)。非参数检验在假设非常不合理(等级有序数)更有效
sapply(UScrime[c("U1","U2")],median)
with(UScrime, wilcox.test(U1,U2,paired = TRUE))
#非参数多于两组的比较
class = state.region
var = state.x77[,c("Illiteracy")]
mydata = as.data.frame(cbind(class,var))
rm(class,var)
install.packages("npmc")
library(npmc)
summary(npmc(mydata),type="BF")
#组件差异的可视化






评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值