# R/es.R In es.dif: Compute Effect Sizes of the Difference

#### Documented in es.ces.des.ees.para.ces.para.des.para.e

```es.d<-function(data1,data2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-inputCheck(data1,data2,alpha,unbiased,vector_out)
if(errorCheck!=0)
{
return(errorCheck)
}
data1<-data1[!is.na(data1)]
data2<-data2[!is.na(data2)]
n1<-length(data1)
n2<-length(data2)
mean1<-mean(data1)
mean2<-mean(data2)
# compute the pooled standard deviation
if(n1!=1)
{
var1<-stats::var(data1)
}
else
{
var1<-0
}
if(n2!=1)
{
var2<-stats::var(data2)
}
else
{
var2<-0
}
return(es.para.d(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out))
}

es.para.d<-function(mean1,mean2,var1,var2,n1,n2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-input.para.Check(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out)
if(errorCheck!=0)
{
return(errorCheck)
}
if(n1+n2<=3)
{
stop("Error:The sum of the sample size must be 4 or larger.")
}
# compute the pooled standard deviation
ss1<-var1*(n1-1)
ss2<-var2*(n2-1)
df<-n1+n2-2
div<-sqrt((ss1+ss2)/df)
# compute the correction factor J
J<-calcJ(df)
# compute d
if(unbiased)
{
d<-(mean1-mean2)/div*J
}
else
{
d<-(mean1-mean2)/div
}
# compute the variance
enye<-n1*n2/(n1+n2)
if(unbiased)
{
para.d<-d
}
else
{
para.d<-d*J
}
var.d<-df/(df-2)*(1/enye+para.d^2)-para.d^2/J^2
if(unbiased)
{
var.d<-var.d*J^2
}
#compute CI
t<-(mean1-mean2)/div*sqrt(enye)
CI<-findCI(alpha,t,df)
if(unbiased)
{
CI<-CI/sqrt(enye)*J
}
else
{
CI<-CI/sqrt(enye)
}
#output
if(vector_out)
{
return(c(d,var.d,CI))
}
else
{
result <- matrix(1, nrow=3,ncol=2)
if(unbiased)
{
result[1,1]<-"Hedges' d:"
}
else
{
result[1,1]<-"Cohen's d:"
}
result[1,2]<-d
result[2,1]<-"variance:"
result[2,2]<-var.d
result[3,1]<-"CI:"
result[3,2]<-paste("[",as.character(CI[1]),",",as.character(CI[2]),"]")
return(result)
}
}

es.e<-function(data1,data2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-inputCheck(data1,data2,alpha,unbiased,vector_out)
if(errorCheck!=0)
{
return(errorCheck)
}
data1<-data1[!is.na(data1)]
data2<-data2[!is.na(data2)]
n1<-length(data1)
n2<-length(data2)
mean1<-mean(data1)
mean2<-mean(data2)
# compute the pooled standard deviation
if(n1!=1)
{
var1<-stats::var(data1)
}
else
{
var1<-0
}
if(n2!=1)
{
var2<-stats::var(data2)
}
else
{
var2<-0
}
return(es.para.e(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out))
}
es.para.e<-function(mean1,mean2,var1,var2,n1,n2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-input.para.Check(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out)
if(errorCheck!=0)
{
return(errorCheck)
}
if((n1==1)||(n2==1))
{
stop("Error: Datum with the sample size 1 cannot be used.")
}
df<-(var1/n1+var2/n2)^2/((var1/n1)^2/(n1-1)+(var2/n2)^2/(n2-1))
div<-sqrt(var1/n1+var2/n2)
enye<-n1*n2/(n1+n2)
# compute the correction factor J
J<-calcJ(df)
# compute e
if(unbiased)
{
e<-(mean1-mean2)/div*J/sqrt(enye)
}
else
{
e<-(mean1-mean2)/div/sqrt(enye)
}
# compute the variance
if(unbiased)
{
para.e<-e
}
else
{
para.e<-e*J
}
var.e<-df/(df-2)*(1/enye+para.e^2)-para.e^2/J^2
if(unbiased)
{
var.e<-var.e*J^2
}
#compute CI
t<-(mean1-mean2)/div
CI<-findCI(alpha,t,df)
if(unbiased)
{
CI<-CI/sqrt(enye)*J
}
else
{
CI<-CI/sqrt(enye)
}
#output
if(vector_out)
{
return(c(e,var.e,CI))
}
else
{
result <- matrix(1, nrow=3,ncol=2)
if(unbiased)
{
result[1,1]<-"Unbiased e:"
}
else
{
result[1,1]<-"Biased e:"
}
result[1,2]<-e
result[2,1]<-"variance:"
result[2,2]<-var.e
result[3,1]<-"CI:"
result[3,2]<-paste("[",as.character(CI[1]),",",as.character(CI[2]),"]")
return(result)
}
}

es.c<-function(data1,data2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-inputCheck(data1,data2,alpha,unbiased,vector_out)
if(errorCheck!=0)
{
return(errorCheck)
}
data1<-data1[!is.na(data1)]
data2<-data2[!is.na(data2)]
n1<-length(data1)
n2<-length(data2)
mean1<-mean(data1)
mean2<-mean(data2)
# compute the standard deviation
if(n1!=1)
{
var1<-stats::var(data1)
}
else
{
var1<-0
}
if(n2!=1)
{
var2<-stats::var(data2)
}
else
{
var2<-0
}
return(es.para.c(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out))
}

es.para.c<-function(mean1,mean2,var1,var2,n1,n2,alpha=0.05,unbiased=TRUE,vector_out=FALSE)
{
errorCheck<-input.para.Check(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out)
if(n1+n2<=3)
{
stop("Error: the sample size of non-constant sample must be 3 or larger.")
}
if(errorCheck!=0)
{
return(errorCheck)
}
# compute the standard deviation
if(n1!=1 && n2==1)
{
ss<-var1*(n1-1)
n<-n1
}
else if(n1==1 && n2!=1)
{
ss<-var2*(n2-1)
n<-n2
}
else if(n1==1 && n2==1)
{
return("Error: both of the data size are 1.")
}
else if(n1!=1 && n2!=1)
{
return("Error: both of the data size are non-1.")
}
df<-n-1
div<-sqrt(ss/df)
# compute the correction factor J
J<-calcJ(df)
# compute c
if(unbiased)
{
c<-(mean1-mean2)/div*J
}
else
{
c<-(mean1-mean2)/div
}
# compute the variance
if(unbiased)
{
para.c<-c
}
else
{
para.c<-c*J
}
var.c<-df/(df-2)*(1/df+para.c^2)-para.c^2/J^2
if(unbiased)
{
var.c<-var.c*J^2
}
#compute CI
t<-(mean1-mean2)/div*sqrt(n-1)
CI<-findCI(alpha,t,df)
if(unbiased)
{
CI<-CI/sqrt(n-1)*J
}
else
{
CI<-CI/sqrt(n-1)
}
#output
if(vector_out)
{
return(c(c,var.c,CI))
}
else
{
result <- matrix(1, nrow=3,ncol=2)
if(unbiased)
{
result[1,1]<-"Unbiased c:"
}
else
{
result[1,1]<-"Biased c:"
}
result[1,2]<-c
result[2,1]<-"variance:"
result[2,2]<-var.c
result[3,1]<-"CI:"
result[3,2]<-paste("[",as.character(CI[1]),",",as.character(CI[2]),"]")
return(result)
}
}

##function to calculate the correction coefficient J
#this function is not for users' direct use
calcJ<-function(m)
{
J<-1
if((m>1)&&(m<=171))
{
J<-gamma(m/2)/sqrt(m/2)/gamma((m-1)/2)
}
else if(m>171)
{
J<-(4*m-4)/(4*m-1)
}
else if(m==1)
{
J<-0
}
else
{
J<-NaN
stop("Unexpected Error: negative degree of freedom.")
}
return(J)
}

##function to search CIs
#this function is not for users' direct use
findCI<-function(alpha,t,df)
{
#check for infinit t
if(is.infinite(t))
{
return(c(NaN,NaN))
}
#find number of digit of alpha/2
x<-alpha/2
digit<-0
while(x %% 1 != 0)
{
x<-x*10
digit<-digit+1
}
#extra digit for more precise computation
digit<-digit+1
#search limits
maxAim<-1-alpha/2
minAim<-alpha/2
#search max
testNcp<--5
nowP<-stats::pt(t,df,testNcp)

while(maxAim>=nowP)
{
testNcp<-testNcp*2
nowP<-stats::pt(t,df,testNcp)
}
minL<-testNcp

#search min
testNcp<-5
while(minAim<=nowP)
{
testNcp<-testNcp*2
nowP<-stats::pt(t,df,testNcp)
}
maxL<-testNcp

#serach CIs

result<-numeric(2)
for(i in 1:2)
{
nowMaxL<-maxL
nowMinL<-minL
nowNcp<-minL+(maxL-minL)/2
nowP<-stats::pt(t,df,nowNcp)
nowP<-signif(nowP,digit)
if(i==1)
{
aimP<-1-alpha/2
}
else
{
aimP<-alpha/2
}
while(nowP!=aimP)
{
if(nowP<aimP)
{
nowMaxL<-nowNcp
}
else
{
nowMinL<-nowNcp
}
nowNcp<-nowMinL+(nowMaxL-nowMinL)/2
nowP<-stats::pt(t,df,nowNcp)
nowP<-signif(nowP,digit)
}
result[i]<-nowNcp
}
return(result)
}

#function to check input
#this function is not for users' direct use
inputCheck<-function(data1,data2,alpha,unbiased,vector_out)
{
if(is.character(data1))
{
return("Error: data1 contains character.")
}
if(is.character(data2))
{
return("Error: data2 contains character.")
}
if(length(data1)==0)
{
return("Error: length of data1 is 0.")
}
if(length(data2)==0)
{
return("Error: length of data2 is 0.")
}
if(is.character(alpha))
{
return("Error: alpha must be a number.")
}
if(alpha<=0||alpha>1)
{
return("Error: alpha must be within (0,1].")
}
if(!is.logical(unbiased))
{
return("Error: unbiased must be TRUE or FALSE.")
}
if(!is.logical(vector_out))
{
return("Error: vector_out must be TRUE or FALSE.")
}
return(0)
}
#function to check input
#this function is not for users' direct use
input.para.Check<-function(mean1,mean2,var1,var2,n1,n2,alpha,unbiased,vector_out)
{
if(is.character(mean1))
{
return("Error: mean1 contains character.")
}
if(is.character(mean2))
{
return("Error: mean2 contains character.")
}
if(is.character(var1))
{
return("Error: var1 contains character.")
}
if(is.character(var2))
{
return("Error: var2 contains character.")
}
if(var1<0)
{
return("Error: var1 must be 0 or larger.")
}
if(var2<0)
{
return("Error: var2 must be 0 or larger.")
}
if(is.character(n1))
{
return("Error: n1 must be a number.")
}
if(is.character(n2))
{
return("Error: n2 must be a number.")
}
if(n1 %% 1!=0||n1<=0)
{
return("Error: n1 must be a natural number.")
}
if(n2 %% 1!=0||n2<=0)
{
return("Error: n2 must be a natural number.")
}
if(is.character(alpha))
{
return("Error: alpha must be a number.")
}
if(alpha<=0||alpha>1)
{
return("Error: alpha must be within (0,1].")
}
if(!is.logical(unbiased))
{
return("Error: unbiased must be TRUE or FALSE.")
}
if(!is.logical(vector_out))
{
return("Error: vector_out must be TRUE or FALSE.")
}
return(0)
}
```

## Try the es.dif package in your browser

Any scripts or data that you put into this service are public.

es.dif documentation built on May 2, 2020, 9:05 a.m.