diff --git a/nimble/model9/fitWildBoarModel9.R b/nimble/model9/fitWildBoarModel9.R
index 80bfb0dceb88f323288e1be74a1bac16e42d273b..7526683c75703203c156e17d58db2bb5ec1011b4 100644
--- a/nimble/model9/fitWildBoarModel9.R
+++ b/nimble/model9/fitWildBoarModel9.R
@@ -1,6 +1,6 @@
# source(here::here("nimble/model9/fitWildBoarModel9.R"))
# rm(list=ls())
-iNSim = 5
+iNSim = 4
source(here::here("src/packages.R"))
source(here::here("src/functions.R"))
diff --git a/nimble/model9/model_definition.R b/nimble/model9/model_definition.R
index 1336f8fa3463c30ad8b316929c56dd37f35b2111..31af07690623e2896aeb741db1f7ba3876b346fb 100644
--- a/nimble/model9/model_definition.R
+++ b/nimble/model9/model_definition.R
@@ -33,14 +33,18 @@ bugsBoarCode = nimbleCode({
logit(scTIntro) ~ dLogitBeta(7, 7) # Introduction date (on 0:1)
logit(pIntroX) ~ dLogitBeta(1, 1) # Introduction location relative to bounding box of whole island
logit(pIntroY) ~ dLogitBeta(1, 1) # Introduction location relative to bounding box of whole island
- logit(pHome) ~ dLogitBeta(2, 2) # Connectivity: pHome=p(wild boar in home cell). (1-pHome)=p(wild boar in neighbouring cell)
+# logit(pHome) ~ dLogitBeta(2, 2) # Connectivity: pHome=p(wild boar in home cell). (1-pHome)=p(wild boar in neighbouring cell)
logit(attractI) ~ dLogitBeta(2, 2) # Relative contact rate, i.e. how attractive is an infectious individual relative to a carcass.
logit(omegaFence) ~ dLogitBeta(2, 2) # Fence efficacy
log(tauP) ~ dnorm(0, sd=3) # Carcass detection rate - passive detection. TESTING THIS PRIOR TO AVOID DISAPPEARING TO ZERO
log(tauA) ~ dLogExp(rate = 1/1e+07) # Carcass detection rate - active detection
log(tauHArmy) ~ dLogExp(rate = 1/1e+07) # Increase in hunting rate around fenced area after day 60
log(tauPArmy) ~ dLogExp(rate = 1/1e+07) # Increase in passive search rate around fenced area after day 60
- log(beta) ~ dLogExp(rate = 1) # The scale of transmission - multiply with the realtive attractivity of Infectious and Carcasses.
+# log(beta) ~ dLogExp(rate = 1) # The scale of transmission - multiply with the realtive attractivity of Infectious and Carcasses.
+ # Changes applied to MCMC4 15h Sunday to attempt reducing beta~pHome redundancy
+ log(beta) ~ dLogExp(rate = 0.1)
+ logit(pHome) ~ dLogitBeta(10, 10)
+ ##
betaI <- beta * attractI
betaC <- beta * (1 - attractI)
tIntro <- tMin - scTIntro * tMin
diff --git a/src/functions.R b/src/functions.R
index 83f7c170756732234b384dc61f5e23fbe5e4da31..be891ab637206dc76d13c21e5d8477869bdb39cf 100644
--- a/src/functions.R
+++ b/src/functions.R
@@ -66,7 +66,7 @@ aggregateWildBoarDataFast <- function(AllBoarData, Hex, stepsPerChunk=1) {
boarData <- AllBoarData %>%
filter(DET!="NAS") %>%
select(!c(DATE.CONF, HOST)) %>% ## DATE.CULLING, ID,
- filter(!is.na(rowHex)) %>% ## DP: there are some NAs around the coast where hexagons and points do not over lap
+ filter(!is.na(rowHex)) %>% ## DP: there are some NAs around the coast where hexagons and points do not over lap
select(DET, DATE.SUSP, rowHex) %>%
st_drop_geometry()
## %>% filter(rowHex <= 3)
@@ -1541,3 +1541,41 @@ firstSampleAboveMean <- function(x) {
}
return(ii)
}
+
+
+scenariosWildBoarReport2 <- function() {
+ flagFence = c(1, 0) ## if 0 -> set logit_omegaFence to -Inf
+ flagHuntZone = c(1, 0) ## if 0 -> set log_tauHArmy to -Inf
+ radiusAS = c(1, 2) ## controls prob of AS in neighbouring cells
+ scenarios = expand.grid(flagFence=flagFence, flagHuntZone=flagHuntZone, radiusAS=radiusAS)
+ return(scenarios)
+}
+
+
+getFarmExposureToInfectiousWildBoar <- function(simIWB_scenario, pigSites, HexAll, Outbreaks, thresh=1) {
+ ## HexAll row index for each pig farm
+ (pigSites$id_HexAll <- as.numeric(st_within(pigSites, HexAll)))
+ ## HexAll row index for each farm
+ (id_HexAllPolygonsWithFarms <- sort(unique(as.numeric(st_within(pigSites, HexAll)))))
+ ## HexAll row index for hexagons with expected cumulative exposure of at least one
+ (id_HexAllCasesInWB <- which(colSums(simIWB_scenario) > thresh))
+ ## Indices of hexagons with exposed farms
+ (id_HexAllWithExposedFarms <- intersect (id_HexAllPolygonsWithFarms, id_HexAllCasesInWB))
+ ## Ranked matrix of most exposed hexagons
+ (HexAllExposed <- t(t(sort(colSums(simIWB_scenario)[id_HexAllWithExposedFarms], decreasing = TRUE))))
+ (HexAllExposed <- cbind(as.integer(sub("hex","",rownames(HexAllExposed))), HexAllExposed))
+ (colnames(HexAllExposed) <- c("idHexAll", "exposure"))
+ HexAllExposed
+ ## Subset of farms and their exposure level
+ (pigSitesAtRisk <- pigSites[is.element(pigSites$id_HexAll, HexAllExposed[,"idHexAll"]),])
+ (pigSitesAtRisk$exposure <- HexAllExposed[ match(pigSitesAtRisk$id_HexAll, HexAllExposed[,"idHexAll"]),"exposure"])
+ (pigSitesAtRisk <- pigSitesAtRisk[order(pigSitesAtRisk$exposure, decreasing=TRUE),])
+ ## identify if they already had cases
+ (id_pig_cases = Outbreaks %>% filter(!is.na(ID)) %>% select(ID) %>% st_drop_geometry())
+ ## Add flag for known cases
+ pigSitesAtRisk$knownCases <- is.element(pigSitesAtRisk$population_id, id_pig_cases$ID)
+ ## Time (day) of peak exposure
+ pigSitesAtRisk$tPeakExposure = apply(simIWB_scenario[, pigSitesAtRisk$id_HexAll], 2, function(x) which(x==max(x))-1)
+ ## Output
+ pigSitesAtRisk
+}
diff --git a/src/plan.R b/src/plan.R
index 76081630b1ef73f76cb3ab231a8b7a259cd1fad1..1108579b236b92a702f10812865303b799ec4e00 100644
--- a/src/plan.R
+++ b/src/plan.R
@@ -730,8 +730,7 @@ plan <- drake_plan(
) %>%
st_drop_geometry() %>%
select(id = population_id, detected, p_pred = p_fomites) %>%
- ## The accuracy of the risk assessment under this threshold is
- ## ridiculous. Better neglect it.
+ ## The accuracy of the risk assessment under this threshold is ridiculous. Better neglect it.
filter(p_pred > 0.001),
.id = "pathway"
) %>%
@@ -783,7 +782,60 @@ plan <- drake_plan(
## WILD BOAR SIMULATION MODEL ####
##%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- ## source(here("nimble/model8/fitWildBoarModel8.R"))
+ ## here("nimble/model1/fitWildBoarModel1.R") ## Estimation for report 1
+ ## here("nimble/simWildBoarFromMCMC.R") ## Simulation for report 1
+
+ ## here("nimble/model9/fitWildBoarModel9.R") ## Estimation for report 2
+ ## here("nimble/model9/simWildBoarFromMCMC9.R") ## Simulation for report 2
+
+ simDirectoryReport2 = here("nimble/model9/MCMC3/sims") # Hopefully we can "upgrade" to MCMC4/sims beofre submission
+ ,
+
+ scenariosWBR2 = scenariosWildBoarReport2()
+ ,
+
+ simIWB_meanI_111 = read.table(
+ file_in(!!here("nimble/model9/MCMC3/sims/Mean_Fence1_Army1_Radius1_nSim5000.csv") ),
+ row.names = 1, col.names = paste0("hex",0:nrow(hexAll))
+ ## !!here("data/TimeSeries_day_80_all.csv"))
+ ),
+
+ ## ## hexAll row index for each pig farm
+ ## (pig_sites$id_hexAll <- as.numeric(st_within(pig_sites, hexAll)))
+ ## ## POINTS of outdoor pig farms
+ ## (pt_outdoorSites <- pig_sites%>%filter(is_outdoor==1))
+ ## ## hexAll row index for outdoor sites
+ ## (id_hexAllPolygonsWithOutdoorFarms <- sort(unique(as.numeric(st_within(pt_outdoorSites, hexAll)))))
+ ## ## hexAll row index for hexagons with expected cumulative exposure of at least one
+ ## (id_hexAllCasesInWB <- which(colSums(simIWB_meanI_111) > 1))
+ ## ## Indices of hexagons with exposed outdoor farms
+ ## (id_hexAllWithExposedFarms <- intersect (id_hexAllPolygonsWithOutdoorFarms, id_hexAllCasesInWB))
+ ## (hexAllExposed <- t(t(sort(colSums(simIWB_meanI_111)[id_hexAllWithExposedFarms], decreasing = TRUE))))
+ ## (hexAllExposed <- cbind(as.integer(sub("hex","",rownames(hexAllExposed))), hexAllExposed))
+ ## (colnames(hexAllExposed) <- c("idHexAll", "exposure"))
+ ## hexAllExposed
+ ## (pigSitesExposure <- pig_sites[is.element(pig_sites$id_hexAll, hexAllExposed[,"idHexAll"]),])
+ ## (pigSitesExposure$exposure <- hexAllExposed[ match(pigSitesExposure$id_hexAll, hexAllExposed[,"idHexAll"]),"exposure"])
+ ## (pigSitesExposure <- pigSitesExposure[order(pigSitesExposure$exposure, decreasing=TRUE),])
+ ## (pigSitesExposure %>% st_drop_geometry)
+ ## ## identify if they already had cases
+ ## (id_pig_cases = outbreaks %>% filter(!is.na(ID)) %>% select(ID) %>% st_drop_geometry())
+ ## ## Add flag for known cases
+ ## pigSitesExposure$knownCases <- is.element(pigSitesExposure$population_id, id_pig_cases$ID)
+ ## ## Time (day) of peak exposure
+ ## pigSitesExposure$tPeakExposure = apply(simIWB_meanI_111[, pigSitesExposure$id_hexAll], 2, function(x) which(x==max(x))-1)
+
+ pigSitesExposure = getFarmExposureToInfectiousWildBoar(simIWB_meanI_111, pig_sites, hexAll, outbreaks, thresh=1)
+ ,
+
+ ## pigSitesExposure %>% filter(is_outdoor==1) %>% select(id_hexAll) %>% st_drop_geometry() %>% unlist() -> idExposedHex
+ ## x <- as.integer(sub("D","",rownames(simIWB_meanI_111)))
+ ## plot(x, simIWB_meanI_111[,1], ylim=c(0,max(simIWB_meanI_111)), typ="n", ylab="E[Infectious Wild Boar]", xlab="Day")
+ ## for (hx in idExposedHex) { # 1:nrow(hexAll)
+ ## lines(x, simIWB_meanI_111[,hx])
+ ## #which(simIWB_meanI_111[,hx]==(simIWB_meanI_111[,hx]))
+ ## }
+
##%%%%%%%%%%%%%%%%%%
## OTHER MODELS ####
diff --git a/src/simulationModel.Rmd b/src/simulationModel.Rmd
index 454a92a6fb3758a7731aa76dee1e95feefcc9935..03bfcc2f5a92ae9fa1d35b357df464a3cbdcbf90 100644
--- a/src/simulationModel.Rmd
+++ b/src/simulationModel.Rmd
@@ -327,10 +327,9 @@ In practice, the weighted average is dominated by a very small number of simulat
Thus the weighted average looks very similar to the most likely simulation -- the main difference being that in the weighted average we see a greater number of pixels
with some (most likely very small) value between zero and one.
-[I can't get drake to embed the video (hosted on Dropbox), but you can see it by clicking here.](https://www.dropbox.com/home/challenge?preview=combined1_iRow5866_-689115.mp4)
-
-[And here is a cmplimentary video showing the weighted standard deviation.](https://www.dropbox.com/s/oboo6f8eur0jxan/standDevI_iRow5866_-689115.mp4?dl=0)
+[I can't get drake to embed the video (hosted on Dropbox), but you can see it by clicking here.](https://www.dropbox.com/home/challenge-ASF?preview=combined1_iRow5866_-689115.mp4)
+[And here is a complimentary video showing the weighted standard deviation.](https://www.dropbox.com/home/challenge-ASF?preview=standDevI_iRow5866_-689115.mp4)
diff --git a/src/submission_2.Rmd b/src/submission_2.Rmd
index 2beeba069530de555445cd11ade904479c829a1a..826cd6ba7579585156fe94f4e0efb2d9d39a7999 100644
--- a/src/submission_2.Rmd
+++ b/src/submission_2.Rmd
@@ -66,8 +66,10 @@ loadd(
pig_sites_pred_D50,
pig_sites_pred_D80,
pig_sites_risks_D80,
+ pigSitesExposure,
tab_farms_pred_D50,
tab_farms_pred_D80,
+ simIWB_meanI_111,
tmAdminPigTypesAndOutbreaks,
tmProportionAgricultureMap,
tmProportionForestCoverMap,
@@ -82,30 +84,30 @@ loadd(
```{r outbreaks_predD50_vs_realised}
-outbreaks_predD50_vs_realised <-
+outbreaks_predD50_vs_realised <-
full_join(
outbreaks_at_D80 %>%
- filter(HOST == "pig herd") %>%
- st_drop_geometry() %>%
+ filter(HOST == "pig herd") %>%
+ st_drop_geometry() %>%
mutate(Infection = ifelse(DATE.SUSP <= 50, "Initial", "New"))
,
- tab_farms_pred_D50 %>%
- select(-size) %>%
- mutate(predicted = p_pred > 0.05) %>%
+ tab_farms_pred_D50 %>%
+ select(-size) %>%
+ mutate(predicted = p_pred > 0.05) %>%
filter(predicted),
by = c(ID = "id")
- ) %>%
+ ) %>%
replace_na(
list(
Infection = "Not infected",
predicted = FALSE
)
- ) %>%
- mutate(Infection = fct_inorder(Infection)) %>%
+ ) %>%
+ mutate(Infection = fct_inorder(Infection)) %>%
left_join(
pig_sites,
by = c(ID = "population_id")
- ) %>%
+ ) %>%
st_sf()
```
@@ -116,11 +118,11 @@ cap <- "Predicted risks (using updated methods) and realised infections in pig s
tm_shape(
admin,
bbox = st_bbox(
- outbreaks_predD50_vs_realised %>%
+ outbreaks_predD50_vs_realised %>%
st_buffer(dist = 5e4)
)
) +
- tm_layout(bg.color = "lightblue") +
+ tm_layout(bg.color = "lightblue") +
# tm_grid(col = "grey", labels.size = 1) +
tm_scale_bar(
breaks = c(0, 50, 100, 150),
@@ -145,9 +147,9 @@ tm_shape(
# tm_dots(col = "FarmType", size=0.3, palette = "viridis") +
tm_shape(
animal_movements %>%
- filter(date >= 50 - detection_delay, date <= 50) %>%
+ filter(date >= 50 - detection_delay, date <= 50) %>%
filter(
- source %in%
+ source %in%
filter(pig_sites_pred_D50, susceptible_wb, !detected)$population_id,
dest %in% farm_pred_inf_ship_D50$dest
)
@@ -161,13 +163,13 @@ tm_shape(
```
```{r unpredicted-sites}
-unpredicted_infections <-
- outbreaks_predD50_vs_realised %>%
+unpredicted_infections <-
+ outbreaks_predD50_vs_realised %>%
filter(!predicted, Infection == "New", is_outdoor == 0)
-unpredicted_infections %>%
- st_drop_geometry() %>%
- select(ID, DET, starts_with("DATE"), FarmType, practice, size) %>%
+unpredicted_infections %>%
+ st_drop_geometry() %>%
+ select(ID, DET, starts_with("DATE"), FarmType, practice, size) %>%
kable(
format = "latex",
booktabs = TRUE,
@@ -179,7 +181,7 @@ unpredicted_infections %>%
```{r susp-contact-case}
-susp_contact_case <-
+susp_contact_case <-
filter(pig_sites_pred_D50, susceptible_wb)[
st_nearest_feature(
unpredicted_infections,
@@ -187,18 +189,18 @@ susp_contact_case <-
),
]
-dist_to_susp_contact_case <-
+dist_to_susp_contact_case <-
st_distance(
unpredicted_infections,
susp_contact_case
) %>% as.numeric
-pig_sites_pred_D80 %>%
+pig_sites_pred_D80 %>%
inner_join(
susp_contact_case %>% st_drop_geometry() %>% select(population_id)
- ) %>%
- st_drop_geometry() %>%
- select(Id = population_id, FarmType, practice, size, "day symptoms" = date_susp) %>%
+ ) %>%
+ st_drop_geometry() %>%
+ select(Id = population_id, FarmType, practice, size, "day symptoms" = date_susp) %>%
kable(
format = "latex",
booktabs = TRUE,
@@ -212,7 +214,7 @@ We did not identify it because it was an __in-door__ site which we had excluded
previous analysis.
However, it is in close proximity (~`r round(dist_to_susp_contact_case, -1)` m) to another infected very small __outdoor__ finishing farm (Table \@ref(tab:susp-contact-case))
-This suggests a fomite transmission pathway, related to proximity to infected
+This suggests a fomite transmission pathway, related to proximity to infected
farms. We include this pathway in the current prediction iteration.
\clearpage
@@ -229,7 +231,7 @@ We consider the risk higher when infected farms are in close proximity and thus
Since infected farms in the viccinity could have not yet been detected, we consider all neighbouring farms weighted by their exposition to infectious wild boar (TODO: probability of infection due to WB?)
-Another change from the previous submission is a correction of the estimation of the probability of infection via animal shipment.
+Another change from the previous submission is a correction of the estimation of the probability of infection via animal shipment.
Previously we considered all the recorded shipments as potential transmission events. Now, we consider only shipments that have taken place in the last 15 days. Earlier shipments can't be transmission events, otherwise the destination farm site would have been detected by now.
@@ -259,20 +261,20 @@ sites_area <- pig_sites_pred_D80 %>% filter(expo_inf_wb > 0) # not only outdoo
```{r}
-distance_to_nearest_infected_site <-
+distance_to_nearest_infected_site <-
st_distance(
- pig_sites_risks_D80 %>%
+ pig_sites_risks_D80 %>%
filter(!detected, susceptible_wb),
- pig_sites_risks_D80 %>%
+ pig_sites_risks_D80 %>%
filter(detected == 1, susceptible_wb),
- ) %>%
+ ) %>%
apply(1, min)
```
```{r risk-fomites, include = FALSE, eval = FALSE}
## Assume a exponentially decreasing exposition kernel
-## with mean of 700 m (the only suspected case) and multiply
+## with mean of 700 m (the only suspected case) and multiply
## by the probability of wild boar infection, taking 1
## for already infected farms.
@@ -282,8 +284,8 @@ distance_to_nearest_infected_site <-
pig_sites_pred_D80 %>%
mutate(
p_fomites = prob_fomites(p_inf, unclass(st_distance(.)))
- ) %>%
- filter(p_fomites > 0.005) %>%
+ ) %>%
+ filter(p_fomites > 0.005) %>%
ggplot(aes(p_inf, p_fomites)) +
geom_point() +
geom_point(
@@ -299,11 +301,11 @@ pig_sites_pred_D80 %>%
# scale_y_log10() +
coord_fixed()
-idx <-
+idx <-
pig_sites_pred_D80 %>%
mutate(
p_fomites = prob_fomites(p_inf, unclass(st_distance(.)))
- ) %>%
+ ) %>%
slice_max(p_fomites, n = 5)
```
@@ -324,11 +326,11 @@ points are outdoor farms. Those already infected and detected by day 80
are marked with a black dot. kde indicates kernel density estimate,
and is used to quantify the exposition level of outdoor farms."
-infected_sites_region <-
+infected_sites_region <-
outbreaks_at_D80 %>%
- filter(HOST == "pig herd") %>%
- select(ID) %>%
- st_set_agr("identity") %>%
+ filter(HOST == "pig herd") %>%
+ select(ID) %>%
+ st_set_agr("identity") %>%
st_intersection(st_union(wb_case_density_D80))
ggplot() +
@@ -396,21 +398,21 @@ ggplot() +
cap <- "Probability of farm infection and detection (curve), estimated from pressence/abscence data (points),
for outdoor farms in the region at risk, as a function of the exposition level."
pig_sites_risks_D80 %>%
- filter(susceptible_wb) %>%
+ filter(susceptible_wb) %>%
mutate(
detected = as.numeric(detected)
) %>%
left_join(
outbreaks_at_D80 %>%
- filter(HOST == "pig herd") %>%
- st_drop_geometry() %>%
+ filter(HOST == "pig herd") %>%
+ st_drop_geometry() %>%
transmute(
population_id = ID,
Infection = ifelse(DATE.SUSP <= 50, "Initial", "New")
),
by = "population_id"
- ) %>%
- replace_na(list(Infection = "Not infected")) %>%
+ ) %>%
+ replace_na(list(Infection = "Not infected")) %>%
ggplot(aes(expo_inf_wb, detected)) +
geom_point(aes(group = population_id, colour = Infection)) +
geom_line(
@@ -458,8 +460,8 @@ plot_area <- st_bbox(
c(
pig_sites_pred_D80 %>% filter(susceptible_wb) %>% st_geometry(),
pig_sites_pred_D80 %>% filter(susceptible_wb) %>% st_geometry() %>% `+`(c(5e4, 0))
- ) %>%
- st_union() %>%
+ ) %>%
+ st_union() %>%
st_buffer(dist = 1e4)
)
@@ -513,7 +515,7 @@ tmAdmin +
tm_shape(
animal_movements %>%
filter(
- source %in%
+ source %in%
filter(pig_sites_pred_D80, susceptible_wb)$population_id,
dest %in% farm_pred_inf_ship_D80$dest)
) +
@@ -531,7 +533,7 @@ tmAdmin +
```{r prob-infection-fomites, fig.cap = cap}
cap <- "Probability of detection due to exposition to fomites during the prediction period (D65 - D95)."
-sites_fomites <-
+sites_fomites <-
pig_sites_pred_D80 %>%
mutate(
p_fomites = prob_fomites(p_inf, unclass(st_distance(.)))
@@ -542,8 +544,8 @@ plot_area <- st_bbox(
c(
sites_fomites %>% st_geometry(),
sites_fomites %>% st_geometry() %>% `+`(c(7e4, 0))
- ) %>%
- st_union() %>%
+ ) %>%
+ st_union() %>%
st_buffer(dist = 1e4)
)
@@ -577,7 +579,7 @@ tm_shape(admin, bbox = plot_area) +
```{r table-farms-pred}
kable(
tab_farms_pred_D80 %>%
- filter(p_pred > 0.05, !detected) %>%
+ filter(p_pred > 0.05, !detected) %>%
select(-detected),
caption = "Probability of infection in undetected farm sites with non-negligible probability (P > 0.05) in the next period.",
format = "latex",
@@ -592,9 +594,9 @@ kable(
```{r prob-infection-all, fig.cap = cap}
cap <- "Probability of infection and detection by any possible route during the prediction period (D65 - D95)."
-sites_all <-
- pig_sites_pred_D80 %>%
- select(-detected, -size, -p_pred) %>%
+sites_all <-
+ pig_sites_pred_D80 %>%
+ select(-detected, -size, -p_pred) %>%
inner_join(
bind_rows(
tab_farms_pred_D80 %>%
@@ -609,8 +611,8 @@ plot_area <- st_bbox(
c(
sites_all %>% st_geometry(),
sites_all %>% st_geometry() %>% `+`(c(7e4, 0))
- ) %>%
- st_union() %>%
+ ) %>%
+ st_union() %>%
st_buffer(dist = 1e4)
)
@@ -649,6 +651,42 @@ tm_shape(admin, bbox = plot_area) +
# Wild boar
+
+```{r table-farm-exposure-IWB-outdoor}
+kable(
+ pigSitesExposure %>% st_drop_geometry() %>% # filter(knownCases == FALSE) %>%
+ filter(is_outdoor==1) %>%
+ select(population_id, size, FarmType, practice, exposure, tPeakExposure, knownCases) ,
+ caption = "Table of outdoor farms most exposed to infectious wild boar, where:
+ exposure_ is the sum, over days 0-140, of the expected number of infectious wild boar within each farm's hexagonal pixel;
+ and tPeakExposure provides the timing (day) with the gretest local density of infectious wild boar.",
+ format = "latex",
+ digits = 3,
+ booktabs = TRUE
+)
+```
+
+
+```{r table-farm-exposure-IWB-indoor}
+kable(
+ pigSitesExposure %>% st_drop_geometry() %>% # filter(knownCases == FALSE) %>%
+ filter(is_outdoor==0) %>%
+ select(population_id, size, FarmType, practice, exposure, tPeakExposure, knownCases)
+ ,
+ caption = "Table of indoor farms located in pixels with greatest exposed to infectious wild boar,
+ where:
+ exposure_ is the sum, over days 0-140, of the expected number of infectious wild boar within each farm's hexagonal pixel;
+ and tPeakExposure provides the timing (day) with the gretest local density of infectious wild boar.",
+ format = "latex",
+ digits = 3,
+ booktabs = TRUE
+)
+```
+
+
+
+
+
\clearpage
# Conclusions and prediction requests
@@ -661,5 +699,5 @@ infected in this period. To a lesser extent ($.1 \leq p \leq .35$), farm sites `
The total number of animals affected in these farms are
`r tab_farms_pred_D80 %>% filter(!detected, p_pred >= .35) %>% pull(size) %>% sum` in the first case and `r tab_farms_pred_D80 %>% filter(!detected, p_pred > .1, p_pred < .35) %>% pull(size) %>% sum` in the second.
-
+
The expected number of affected animals (in addition to those from already detected farms) is `r (exp_nanimals = tab_farms_pred_D80 %>% filter(!detected) %>% with(., sum(p_pred * size))) %>% round()` out of `r round((total_animals = pig_sites_risks_D80 %>% filter(!detected) %>% pull(size) %>% sum) %>% "/"(1e6))` M heads (`r round(exp_nanimals / total_animals * 100, 2)`%).