prop.trend.test.moulton: prop.trend.test.moulton

In tbgitoo/moultonTools: moultonTools

Description

Proportion trend test, with Moulton correction for clustering.

Usage

prop.trend.test.moulton(binaryResponse, score, cluster=seq_along(binaryResponse),... )

Arguments

binaryResponse	A vector of TRUE and FALSE values ("logical")
score	A numerical vector of the same length as binaryResponse. This variable is used to evaluate the proportion of succesful (TRUE) outcomes for each value of score; these proportions are the basis for prop.trend.test invoked by this function.
cluster	A vector indicating the cluster identities, of the same length as binaryResponse
...	Additional arguments to be passed to the underlying moulton_factor function, particularly a method argument to use different methods for the estimation of intraclass correlation. See also ICC

.

Details

The function evaluates the proportions from the binary response variable (binaryResponse) categorized using the score variable. It then invokes prop.trend.test with these proportions. Then, the Moulton correction is applied, by means of a Moulton factor calculated by moulton_factor based on the cluster column

Value

An object of class "htest", see prop.trend.test. In addition, the Moulton factor used is indicated in a "moulton_factor" field.

Author(s)

Marina + Thomas Braschler

Examples

# Example 1: Construction of clustered variables and testing


# The primary question about Moulton correction is false positive detection

## Not run: 

one_round_simulation<-function(p0=0.3,sd_cluster=1,n_p=3,n_cluster_per_p=2,n_per_cluster=10)
{


x0 = log(p0/(1+p0))  # To scale the interval 0,1 to -inf to +inf


xtot = vector(mode="numeric",length=0)
score =  vector(mode="numeric",length=0)
cluster_id =  vector(mode="numeric",length=0)


for(ind_score in 1:n_p)
{


score = c(score, rep(ind_score,n_per_cluster*n_cluster_per_p))

	for(ind_cluster in 1:n_cluster_per_p)
	{
		cluster_id = c(cluster_id, rep((ind_score-1)*n_cluster_per_p+ind_cluster,n_per_cluster))
		xtot = c(xtot,rep(x0+rnorm(1,sd=sd_cluster),n_per_cluster))
	}

}

ptot=exp(xtot)/(1+exp(xtot))


response = vector(mode="numeric",length=0)


for(ind in 1:length(ptot))
{
	response=c(response,rbinom(1,1,ptot[ind]))
}


x_prop_trend_test = aggregate(response ~ score, FUN=sum)$response
n_prop_trend_test = aggregate(response ~ score, FUN=length)$response

p_vector=vector(mode="numeric",length=2)

names(p_vector)=c("no.correction","moulton")

# There are occasionally some occurrences where the proportions generate warnings, suppress it here. This happens if a perfect with is possible
# for example if x_prop_trend_test is constant or if makes a regularly ascending or descending series as in (2,5,8) out of (20,20,20). This warning here
# is irrelevant, suppress
p_vector["no.correction"]=suppressWarnings(prop.trend.test(x_prop_trend_test,n_prop_trend_test)$p.value)
p_vector["moulton"]=suppressWarnings(prop.trend.test.moulton(response, score, cluster_id )$p.value)

return(p_vector)


}

p_matrix=matrix(nrow=500,ncol=2)

colnames(p_matrix)=c("no.correction","moulton")


for(ind in 1:(dim(p_matrix)[1]))
{

	p_matrix[ind,]=one_round_simulation(p0=0.3,sd_cluster=0.5,n_p=
	3,n_cluster_per_p=2,n_per_cluster=40)


}



cat(paste("Fraction of false positives at alpha = 0.05 WITHOUT Moulton correction: ",sum(p_matrix[,1]<0.05)/(dim(p_matrix)[1]),sep=""))

cat(paste("Fraction of false positives at alpha = 0.05 WITH Moulton correction: ",sum(p_matrix[,2]<0.05)/(dim(p_matrix)[1]),sep=""))


## End(Not run)

tbgitoo/moultonTools documentation built on Nov. 8, 2021, 6:15 p.m.