library(groupTesting)
## To illustrate 'prop.gt', we use data simulated by
## the R functions 'hier.gt.simulation' and 'array.gt.simulation'.
## The simulated data-structures are consistent
## with the data-structure required for 'gtData'.
## Example 1: MLE from 3-stage hierarchical group testing data.
## The data used is simulated by 'hier.gt.simulation'.
N <- 90 # Sample size
S <- 3 # 3-stage hierarchical testing
psz <- c(6,2,1) # Pool sizes used in stages 1-3
Se <- c(.95,.95,.98) # Sensitivities in stages 1-3
Sp <- c(.95,.98,.96) # Specificities in stages 1-3
assayID <- c(1,2,3) # Assays used in stages 1-3
p.t <- 0.05 # The TRUE parameter to be estimated
# Simulating data:
set.seed(123)
gtOut <- hier.gt.simulation(N,p.t,S,psz,Se,Sp,assayID)$gtData
# Running the EM algorithm:
pStart <- p.t + 0.2 # Initial value
res <- prop.gt(p0=pStart,gtData=gtOut,covariance=TRUE,
nburn=2000,ngit=5000,maxit=200,tol=1e-03,
tracing=TRUE,conf.level=0.95)
# Estimation results:
# > res
# $param
# [1] 0.05158
# $covariance
# [,1]
# [1,] 0.0006374296
# $iterUsed
# [1] 4
# $convergence
# [1] 0
# $summary
# Estimate StdErr 95%lower 95%upper
# prop 0.052 0.025 0.002 0.101
## Example 2: MLE from two-dimensional array testing data.
## The data used is simulated by 'array.gt.simulation'.
N <- 100 # Sample size
protocol <- "A2" # 2-stage array without testing the initial master pool
n <- 5 # Row/column size
Se <- c(0.95, 0.95) # Sensitivities
Sp <- c(0.98, 0.98) # Specificities
assayID <- c(1, 1) # The same assay in both stages
p.true <- 0.05 # The TRUE parameter to be estimated
# Simulating data:
set.seed(123)
gtOut <- array.gt.simulation(N,p.true,protocol,n,Se,Sp,assayID)$gtData
# Fitting the model:
pStart <- p.true + 0.2 # Initial value
res <- prop.gt(p0=pStart,gtData=gtOut,covariance=TRUE)
print(res)
# \donttest{
## Example 3: MLE from non-overlapping initial pooled responses.
## The data used is simulated by 'hier.gt.simulation'.
## Note: With initial pooled responses, our MLE is equivalent
## to the MLE in Litvak et al. (1994) and Liu et al. (2012).
N <- 1000 # Sample size
psz <- 5 # Pool size
S <- 1 # 1-stage testing
Se <- 0.95 # Sensitivity
Sp <- 0.99 # Specificity
assayID <- 1 # Assay used for all pools
p.true <- 0.05 # True parameter
set.seed(123)
gtOut <- hier.gt.simulation(N,p.true,S,psz,Se,Sp,assayID)$gtData
pStart <- p.true + 0.2 # Initial value
res <- prop.gt(p0=pStart,gtData=gtOut,
covariance=TRUE,nburn=2000,ngit=5000,
maxit=200,tol=1e-03,tracing=TRUE)
print(res)
## Example 4: MLE from individual (one-by-one) testing data.
## The data used is simulated by 'hier.gt.simulation'.
N <- 1000 # Sample size
psz <- 1 # Pool size 1 (i.e., individual testing)
S <- 1 # 1-stage testing
Se <- 0.95 # Sensitivity
Sp <- 0.99 # Specificity
assayID <- 1 # Assay used for all pools
p.true <- 0.05 # True parameter
set.seed(123)
gtOut <- hier.gt.simulation(N,p.true,S,psz,Se,Sp,assayID)$gtData
pStart <- p.true + 0.2 # Initial value
res <- prop.gt(p0=pStart,gtData=gtOut,
covariance=TRUE,nburn=2000,
ngit=5000,maxit=200,
tol=1e-03,tracing=TRUE)
print(res)
## Example 5: Using pooled testing data.
# Pooled test outcomes:
Z <- c(1, 0, 1, 0, 1, 0, 1, 0, 0)
# Pool sizes used:
psz <- c(6, 6, 2, 2, 2, 1, 1, 1, 1)
# Pool-specific Se & Sp:
Se <- c(.90, .90, .95, .95, .95, .92, .92, .92, .92)
Sp <- c(.92, .92, .96, .96, .96, .90, .90, .90, .90)
# Assays used:
Assay <- c(1, 1, 2, 2, 2, 3, 3, 3, 3)
# Pool members:
Memb <- rbind(
c(1, 2, 3, 4, 5, 6),
c(7, 8, 9, 10, 11, 12),
c(1, 2, -9, -9, -9, -9),
c(3, 4, -9, -9, -9, -9),
c(5, 6, -9, -9, -9, -9),
c(1,-9, -9, -9, -9, -9),
c(2,-9, -9, -9, -9, -9),
c(5,-9, -9, -9, -9, -9),
c(6,-9, -9, -9, -9, -9)
)
# The data-structure suited for 'gtData':
gtOut <- cbind(Z, psz, Se, Sp, Assay, Memb)
# Fitting the model:
pStart <- 0.10
res <- prop.gt(p0=pStart,gtData=gtOut,
covariance=TRUE,nburn=2000,
ngit=5000,maxit=200,
tol=1e-03,tracing=TRUE)
print(res)
# }
Run the code above in your browser using DataLab