For a sampling scheme designed to substantiate freedom from disease the function computes the a-posteriori alpha-error, i.e., the actual alpha-error based on the drawn sample, when the population is stratified into risk groups for which the relative infection risk is known.

1 2 3 | ```
computeAposterioriErrorRiskGroups(alphaErrorVector,
groupVec, groupLevels, nPopulationVec, nRelRiskVec,
nDiseased, method = "default")
``` |

`alphaErrorVector` |
Numeric vector. Alpha-error (between 0 and 1) of each herd in the sample. |

`groupVec` |
Character vector. Risk group to which each of the
herds in the sample belongs. Must have the same length
(and order) as |

`groupLevels` |
Character vector. (Unique) levels/names of the risk group. Defines the order
of the values in |

`nPopulationVec` |
Integer vector. Population sizes of the risk groups. Must have the same length
(and order) as |

`nRelRiskVec` |
Numeric vector. (Relative) infection risks of the risk groups. Must have
the same length (and order) as |

`nDiseased` |
Integer. Number of diseased herds in the population according to the design prevalence. |

`method` |
Character string. "exact" for exact error, "approx" for approximation (recommended for nDiseased > 7). |

The exact evaluation of the alpha-error is computationally complex, due
to combinatirical issues. In order to increase effectivity parts of the code
were implemented in C. Still, for `nDiseased`

> 3 the computation may take
very long and it is generally not recommended to use the exact calculation. Rather
the approximation should be used for `nDiseased`

> 3.

The return value is the a-posteriori alpha-error based on the sample at hand (numeric scalar).

Ian Kopacka <ian.kopacka@ages.at>

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 | ```
data(sheepData)
sheepData$size <- ifelse(sheepData$nSheep < 30, "small", "large")
riskValueData <- data.frame(riskGroup = c("small", "large"),
riskValues = c(1,2))
mySurvey <- surveyData(nAnimalVec = sheepData$nSheep,
riskGroupVec = sheepData$size,
riskValueData = riskValueData,
populationData = sheepData, designPrevalence = 0.002,
alpha = 0.05, intraHerdPrevalence = 0.13,
diagSensitivity = 0.9, costHerd = 30, costAnimal = 7.1)
## Limited sampling with risk groups:
myLtdSamplingRG <- ltdSampling(survey.Data = mySurvey, sampleSizeLtd = 20,
nSampleFixVec = NULL, probVec = c(1,2))
## Draw sample manually:
set.seed(20110708)
x <- myLtdSamplingRG
indexSampleRG <- sapply(seq(along = x@surveyData@riskValueData$riskGroup),
function(ii){
riskGroup <- as.character(x@surveyData@riskValueData$riskGroup[ii])
nSampleRG <- x@nHerdsPerRiskGroup[riskGroup]
indexRG <- which(x@surveyData@riskGroupVec == riskGroup)
indexOut <- sample(x = indexRG, size = nSampleRG, replace = FALSE)
})
indexSample <- sort(Reduce(function(x,y) c(x,y), indexSampleRG))
## Compute the a-posteriori alpha error:
alphaErrorVector <- computeAlpha(nAnimalVec = x@surveyData@nAnimalVec[indexSample],
method = "limited", sampleSizeLtd = x@sampleSizeLtd,
intraHerdPrevalence = x@surveyData@intraHerdPrevalence,
diagSensitivity = x@surveyData@diagSensitivity, diagSpecificity = 1)
## Determine the number of herds in each risk group:
riskValueDf <- mySurvey@riskValueData[,1:2]
names(riskValueDf) <- c("riskGroup", "riskValues")
riskValueDf$riskGroup <- as.character(riskValueDf$riskGroup)
riskValueDf$id <- seq(along = riskValueDf[,1])
riskGroupTab <- table(mySurvey@riskGroupVec)
riskGroupDf <- data.frame(riskGroup = as.character(names(riskGroupTab)),
nPopulation = as.vector(riskGroupTab))
riskValueDf <- merge(x = riskValueDf, y = riskGroupDf, by = "riskGroup",
sort = FALSE)
riskValueDf <- riskValueDf[order(riskValueDf$id, decreasing = FALSE),]
aPostAlphaRG <- computeAposterioriErrorRiskGroups(alphaErrorVector = alphaErrorVector,
groupVec = x@surveyData@riskGroupVec[indexSample],
groupLevels = riskValueDf$riskGroup,
nPopulationVec = riskValueDf$nPopulation,
nRelRiskVec = riskValueDf$riskValues,
nDiseased = max(round(length(x@surveyData@nAnimalVec)*x@surveyData@designPrevalence),1),
method = "approx")
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.