Report of ABC for Irish hidden burden model

1 Aim

To trial the model and ABC using target metrics from real data.

2 Methods

2.1 Population

Sample of 1500 breakdowns from all herds (including GIF tested) with breakdowns ending between 2015 and 2018 and with at least three years of followup. Max herd size allowed was 600 and trader herds were excluded.

2.2 Model

Allows density dependence (UK herd size midpoint of 165 was used, true midpoint would be ~ 300). Full model of GIF testing and first case factor (tester behaviour) included.

2.3 ABC

Target metrics are:

• % breakdowns with GIF testing (4 or more case threshold and probability based on regression output);
• % breakdowns with skin test reactors after first test of breakdown;
• % first breakdowns detected by slaughterhouse surveillance;
• % breakdowns with recurrence within 1 year, 2 years and 3 years;
• Binned breakdown duration;
• Binned count of skin reactors in the first test of breakdown;
• Binned count of total skin reactors throughout the breakdown;
• Binned count of GIF positive reactors.

Particles were rejected completely if they did not fulfill the criteria below:

• too few breakdowns;
• breakdowns running beyond 4 years;
• too little recurrence;
• too much recurrence;
• too few slaughterhouse initiated breakdowns;
• too few GIF tested breakdowns;
• too many breakdowns with greater than 20 cases;
• too many breakdowns with greater than 5 slaughterhouse cases.

The 75th percentile of Round 1a likelihoods was used as the threshold for Round 1.However, this ran more easily than I expected and I had more Round 1 particles than I needed, so I restricted to about the 56th percentile of my Round1 KL distances (leaving me with 516 Particles for Round 2). Subsequent round thresholds were 95% of previous KLs, and I made the “discard particles” rules slightly more stringent as rounds progressed.

3 Results

3.1 Posterior predictions

3.1.1 Summary statistics

3.1.1.1 Breakdown duration

``` ```{=html}

(\#tab:four)Summary of breakdown duration.

status	min	q25	median	q75	max	mean	sd
data	54	134	147	174.00	1,293	170.4387	79.54322
low_kl	0	133	149	203.25	1,034	181.1300	93.72183
model	0	134	151	213.00	1,443	187.6368	102.17249

Figure 3.1: First breakdown duration in days.

3.1.2 Proportion of breakdowns with greater than one case

Figure 3.2: Proportion of breakdowns with recurrence within two years, stratified by herd size category.

3.1.3 Breakdown recurrence within three years

Figure 3.3: Proportion of breakdowns with recurrence within two years, stratified by herd size category.

3.1.4 Proportion of breakdowns initiated by lesion detection at routine slaughter

Figure 3.4: Proportion of breakdowns initiated by slaughterhouse detection, stratified by herd size category.

3.2 Averaged outputs for posterior predictions stratified by herd size

For the 100 particle sample from final round, for each herd ID, the 100 outputs were averaged into single outputs.

The black dots are the data and the red is model output.

3.3 Posterior predictive plots stratified by herd management system

Management system may give some hint about drivers of transmission and detection. For example, a subset of dairy animals may have longer residency times. Fattener herds will have a larger in degree and will have more moves to slaughter.

3.3.1 Breakdown duration

3.3.1.1 Median and IQR breakdown duration

Figure 3.5: Median and IQR of first breakdown duration in days, stratified by herd management practice.

3.3.1.2 Proportion breakdowns of greater than 240 days in duration

Figure 3.6: Proportion of breakdowns greater than 240 days in duration. Please note the confidence intervals for the simulations are atrificially small as these are based on output from 100 repetitions of the model with the true fixed or post - ABC median parameters.

3.3.2 Count of cases in breakdown

Figure 3.7: Median and interquartile ranges for count of cases during a breakdown. The black line is the median and the boxes encompass the IQR.

3.3.3 Proportion of breakdowns with greater than one case

Figure 3.8: Proportion of breakdowns with recurrence within two years, stratified by herd management practice.

3.3.4 Breakdown recurrence within three years

Figure 3.9: Proportion of breakdowns with recurrence within two years, stratified by herd management practice.

3.3.5 Proportion of breakdowns initiated by lesion detection at routine slaughter

Figure 3.10: Proportion of breakdowns initiated by slaughterhouse detection, stratified by herd management practice.

3.4 Parameters

3.4.1 Individual parameter plots

Round1a is not quite the prior anymore, as particles not fulfilling certain criteria get discarded.

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

3.5 ABC characteristics

3.5.1 Sample weights in last round

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Figure 3.11: Weights of Particles from last round available.

3.5.2 KL distances in Round1a

## No id variables; using all as measure variables

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

3.5.3 KL distances in last Round

## No id variables; using all as measure variables

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

3.5.4 Particle with lowest KL distance versus weighted sample of 100 particles

``` ```{=html}

(\#tab:three)Parameter comparisons.

Parameter	low KL sample	samples
beta	1.40415752228514e-05	1.23684135691347e-05
skin_tweak	0.468165076566237	1.74094361126135
first_case_factor	10.8021136336427	12.3157518766988
gif_cutoff	-0.843674132443083	-0.248546090929861
chi	0.310870645978277	0.666148513586836
move_chi	0.00469916764200435	0.0317773959996895
sl_sens	0.821010869030637	0.50194885852158
test_occult	45.0359536019094	43.6704673201604
dens_q_param	0.844756303197237	0.867217091896389
sev_cutoff	2.46816507656624	3.74094361126135
std_cutoff	4.46816507656624	5.74094361126135
sev_fc	13.270278710209	16.0566954879601
std_fc	15.270278710209	18.0566954879601
gif_sens	0.92522	0.90829
sev_sens	0.80741	0.75135
std_sens	0.71633	0.65442
sev_sens_fc	0.3059	0.21377
std_sens_fc	0.23696	0.16273
gif_sp	0.82325	0.88349
sev_sp	0.99404	0.99928
std_sp	0.99985	0.99999
sev_sp_fc	1	1
std_sp_fc	1	1

3.5.5 Effective sample size

list_ess <- numeric(length = length_rounds)

for (i in 1:length_rounds)
{
  list_ess[i] <-  ESS(list_rounds[[i]][[3]])
}

plot(list_ess, type = "l", ylim = c(0, Particles))

3.6 Priors

Where possible, uninformative priors were used in this ABC.

``` ```{=html}

(\#tab:one)Prior limits.

parameter	lower_prior	upper_prior
beta	0.000	0.0001818182
skin_tweak	-4.000	8.0000000000
first_case_factor	0.000	15.0000000000
gif_cutoff	-6.000	4.0000000000
chi	0.001	5.0000000000
move_chi	0.001	1.0000000000
sl_sens	0.000	1.0000000000
test_occult	14.000	84.0000000000
dens_q_param	0.000	1.0000000000

## No id variables; using all as measure variables

## $`1`

## 
## $`2`

## 
## $`3`

## 
## $`4`

## 
## attr(,"class")
## [1] "list"      "ggarrange"