《R语言与统计分析》-参数的假设检验

####方差已知时，求单正态总体均值置信区间-z.test####
####单正态总体均值假设检验-z.test####
#R无内置函数，自己写#
z.test<-function(x,n,sigma,alpha,u0=0,alternative="two.sided"){
  options(digits = 4)
  result<-list()
  mean<-mean(x)
  z<-(mean-u0)/(sigma/sqrt(n))
  p<-pnorm(z,lower.tail = FALSE)
  result$mean<-mean
  result$z<-z
  result$p.value<-p
  if(alternative=="two.sided")
    result$p.value<-2*pnorm(abs(z),lower.tail = FALSE)
  else if(alternative=="greater")
    result$p.value<-pnorm(z)
  result$conf.int<-c(
    mean-sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n),
    mean+sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n)
  )
  result
}

####方差未知时，求单正态总体均值置信区间-t.test####
t.test(x, y = NULL,
       alternative = c("two.sided", "less", "greater"),
       mu = 0, paired = FALSE, var.equal = FALSE,
       conf.level = 0.95, ...)
#单样本时，y=NULL，单样本t检验
#二样本时，x，y，二样本t检验
#alternative = c(上下界, 置信上限, 置信下限)

####求单正态总体方差的置信区间-chisq.var.test####
#R无内置函数，自己写#
chisq.var.test<-function(x,var,alpha,alternative="two.sided"){
  options(digits = 4)
  result<-list()
  n<-length(x)
  v=var(x)
  result$var<-v
  chi2<-(n-1)*v/var
  result$chi2<-chi2
  p<-pchisq(chi2,n-1)
  result$p.value<-p
  if(alternative=="less")
    result$p.value<-pchisq(chi2,n-1,lower.tail = F)
  else if(alternative=="two.sided")
    result$p.value<-2*min(pchisq(chi2,n-1),
                          pchisq(chi2,n-1,lower.tail = F))
  result$conf.int<-c(
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = F),
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = T)
  )
  result
}

####两总体均值差的区间估计####
####两方差已知-two.sample.ci####
#R无内置函数，自己写#
two.sample.ci<-function(x,y,conf.level=0.95,sigma1,sigma2){
  options(digits = 4)
  m=length(x);n=length(y)
  xbar=mean(x)-mean(y)
  alpha=1-conf.level
  zstar=qnorm(1-alpha/2)*(sigma1/m+sigma2/n)^0.5
  xbar+c(-zstar,+zstar)
}

####量方差比的置信区间-var.test####
var.test(x, y, ratio = 1,
         alternative = c("two.sided", "less", "greater"),
         conf.level = 0.95, ...)

####单总体比率p的区间估计-prop.test/binom.test####
prop.test(x, n, p = NULL,
          alternative = c("two.sided", "less", "greater"),
          conf.level = 0.95, correct = TRUE)
#correct为是否连续性校正

binom.test(x, n, p = 0.5,
           alternative = c("two.sided", "less", "greater"),
           conf.level = 0.95)
#应用于二项分布

####两总体比率差的区间估计-prop.test/ratio.ic####
smokers  <- c( 83, 90, 129, 70 )
patients <- c( 86, 93, 136, 82 )
prop.test(smokers, patients)

#可自己写#
ratio.ic<-function(x,y,n1,n2,conf.level=0.95){
  xbar1=x/n1;xbar2=y/n2
  xbar=xbar1-xbar2
  alpha=1-conf.level
  zstar=qnorm(1-alpha/2)*(xbar1*(1-xbar1)/n1+xbar2*(1-xbar2)/n2)^0.5
  xbar+c(-zstar,+zstar)
}

####估计正态总体均值时样本容量的确定####
####总体方差已知####
#R无内置函数，自己写#
size.norm1 <- function(d,var,conf.level=0.95) {
  alpha=1-conf.level
  ((qnorm(1-alpha/2)*var^0.5)/d)^2
}
#d是可允许的误差

####总体方差未知####
#R无内置函数，自己写#
size.norm2<-function(s,alpha,d,m){
  t0<-qt(alpha/2,m,lower.tail = F)
  n0<-(t0*s/d)^2
  t1<-qt(alpha/2,n0,lower.tail = F)
  n1<-(t1*s/d)^2
  while(abs(n1-n0)>0.5){
    n0<-(qt(alpha/2,n1,lower.tail = F)*s/d)^2
    n1<-(qt(alpha/2,n0,lower.tail = F)*s/d)^2
  }
  n1
}
#m是一个事先给定的很大的数字
#s是样本标准差
#d是可允许的误差

####估计比率p时的样本容量确定####
#R无内置函数，自己写#
size.bin<-function(d,p,conf.level=0.95){
  alpha=1-conf.level
  ((qnorm(1-alpha/2))/d)^2*p*(1-p)
}

####正态总体均值的假设检验####
####方差已知时-z.test####
z.test<-function(x,n,sigma,alpha,u0=0,alternative="two.sided"){
  options(digits = 4)
  result<-list()
  mean<-mean(x)
  z<-(mean-u0)/(sigma/sqrt(n))
  p<-pnorm(z,lower.tail = FALSE)
  result$mean<-mean
  result$z<-z
  result$p.value<-p
  if(alternative=="two.sided")
    result$p.value<-2*pnorm(abs(z),lower.tail = FALSE)
  else if(alternative=="greater")
    result$p.value<-pnorm(z)
  result$conf.int<-c(
    mean-sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n),
    mean+sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n)
  )
  result
}

####方差未知时-t.test####
t.test(x, y = NULL,
       alternative = c("two.sided", "less", "greater"),
       mu = 0, paired = FALSE, var.equal = FALSE,
       conf.level = 0.95, ...)

####方差检验-chisq.var.test####
#R无内置函数，自己写#
chisq.var.test<-function(x,var,alpha,alternative="two.sided"){
  options(digits = 4)
  result<-list()
  n<-length(x)
  v=var(x)
  result$var<-v
  chi2<-(n-1)*v/var
  result$chi2<-chi2
  p<-pchisq(chi2,n-1)
  result$p.value<-p
  if(alternative=="less")
    result$p.value<-pchisq(chi2,n-1,lower.tail = F)
  else if(alternative=="two.sided")
    result$p.value<-2*min(pchisq(chi2,n-1),
                          pchisq(chi2,n-1,lower.tail = F))
  result$conf.int<-c(
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = F),
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = T)
  )
  result
}

####两正态总体的参数检验####
####均值的比较-t.test####
t.test(x, y = NULL,
       alternative = c("two.sided", "less", "greater"),
       mu = 0, paired = FALSE, var.equal = FALSE,
       conf.level = 0.95, ...)

####方差的比较-F检验-var.test####
var.test(x, y, ratio = 1,
         alternative = c("two.sided", "less", "greater"),
         conf.level = 0.95, ...)

####成对数据检验-t.test/onesamp####
t.test(x, y = NULL,
       alternative = c("two.sided", "less", "greater"),
       mu = 0, paired = T, var.equal = FALSE,
       conf.level = 0.95, ...)
#paired=T代表成对数据

library(DAAG)
onesamp(dset, x="unsprayed", y="sprayed", xlab=NULL, ylab=NULL,
        dubious=NULL, conv=NULL, dig=2, ...)
#dset为两列的矩阵、数据框；x为处于predictor的列名；y为处于response的列名

####单样本比率检验####
####比率p的精确检验-F检验-binom.test####
binom.test(x, n, p = 0.5,
           alternative = c("two.sided", "less", "greater"),
           conf.level = 0.95)

####比率p的近似检验(大样本,n>30)-prop.test####
prop.test(x, n, p = NULL,
          alternative = c("two.sided", "less", "greater"),
          conf.level = 0.95, correct = TRUE)

####两样本比率检验-prop.test####
prop.test(x, n, p = NULL,
          alternative = c("two.sided", "less", "greater"),
          conf.level = 0.95, correct = TRUE)

smokers  <- c( 83, 90, 129, 70 )
patients <- c( 86, 93, 136, 82 )
prop.test(smokers, patients)