Sampling Plan 3D.Rmd

---
title: 'Sampling: Sampling Plan 3D'
author: "Xianbin Cheng"
date: "February 12, 2019"
output:
  html_document: default
---

# Objective

  * Add a 3D component to the sampling plan module

# Method

1. Load libraries and source R code.

```{r, warning = FALSE, message = FALSE}
source(file = "Sampling_libraries.R")
source(file = "Sampling_contamination.R")
source(file = "Sampling_visualization.R")
source(file = "Sampling_plan.R")
source(file = "Sampling_plan_3d.R")
library(plotly)
```

```{r}
sessionInfo()
```

2. List important parameters from `Sampling_contamination.R`.

* `n_contam` = the number of contamination points  
* `x_lim` = the limits of the x-axis  
* `y_lim` = the limits of the y-axis  
* `z_lim` = the limits of the z-axis
* `x` = the horizontal coordinate of the contamination center, which follows a uniform distribution by default
* `y` = the vertical coordinate of the contamination center, which follows a uniform distribution by default
* `z` = the vertical coordinate of the contamination center, which follows a uniform distribution by default
* `cont_level` = a vector that indicates the mean contamination level (logCFU/g or logCFU/mL) and the standard deviation in a log scale, assuming contamination level follows a log normal distribution $ln(cont\_level)$~$N(\mu, \sigma^2)$. 

** Mode 1: Discrete Spread** 

* `n_affected` = the number of affected kernels near the contamination spot, which follows a Poisson distribution ($Pois(\lambda = 5)$)   
* `covar_mat` = covariance matrix of `x` and `y`, which defines the spread of contamination. Assume the spread follows a 2D normal distribution with var(X) = 0.25, var(Y) = 0.25 and cov(X, Y) = 0   

** Model 2: Continuous Spread**

We do not consider such type of spread in a 3D space. 

```{r}
## The input parameters
n_contam = 3
x_lim = c(0, 30)
y_lim = c(0, 10)
z_lim = c(0, 2.5)
lims = list(xlim = x_lim, ylim = y_lim, zlim = z_lim)
cont_level = c(7, 1)
spread = "discrete"

### Discrete
n_affected = 2
covar_mat = make_covar_mat(spread = spread, varx = 0.25, vary = 0.25, varz = 0.25, covxy = 0, covxz = 0, covyz = 0)

### Continuous
spread_radius = 1
LOC = 10^(-3)
fun = "exp"
```

```{r}
# Run this if we want reproducibility
set.seed(123)
```

```{r}
# Generate the coordinates of contamination points
contam_xy = sim_contam_new(n_contam = n_contam, lims = lims, spread = spread, covar = covar_mat, n_affected = n_affected, spread_radius = spread_radius, cont_level = cont_level) 
```

```{r}
## Basic info of the contamination simulation
str(contam_xy)

summary(contam_xy$label)
```

3. Define parameters for generating sampling plans.

* `n_sp` = the number of sampling points
* `sp_radius` = the radius (m) of a circular region around the sample point. (Only applicable to **Mode 1: Discrete Spread**)
* `n_strata` = the number of strata (applicable to *2D and 3D Stratified random sampling*)
* `by` = the side along which the field is divided into strata. It is either "row" or "column" (applicable to *2D Stratified random sampling*) **OR** the side along which a sample is taken every k steps (applicable to *2D Systematic sampling*).
* `method_sp` = sampling strategy, including `srs`, `strs`, `ss`.

**Mode 1: Discrete Spread**   
* `m_kbar` = averaged kernel weight (g)
* `m_sp` = the analytical sample weight (25 g)
* `conc_good` = concentration of toxin in healthy kernels

```{r}
n_sp = 15
sp_radius = 0.5
n_strata = 5
by = "row"
method_sp = "ss"

# m_kbar = 0.3
# m_sp = 25
# conc_good = 0.7 # We could use gen_elisa_af() from STD_simulation.R
```

4. Generate a sampling plan.

* Systematic sampling:

  1) Truck: 
  
    * `depth_ft` = the depth of corn in feet. It could be >= 4 ft or < 4 ft, each corresponding to a different sampling pattern. By default, ft2m(depth_ft) == zlim[2]
  
  2) Barge
  
  3) Hopper car:
  
    * 3-compartment trough type
    
    * 2-compartment open top type


```{r}
## Systematic sampling
# truck
sp_xy1 = sim_plan_new(method_sp = method_sp, spread = spread, lims = append(lims[1:2], values = list(zlim = c(0,ft2m(3.5)))), radius = sp_radius, container = "truck")
sp_xy2 = sim_plan_new(method_sp = method_sp, spread = spread, lims = append(lims[1:2], values = list(zlim = c(0,ft2m(5)))), radius = sp_radius, container = "truck")

# barge
sp_xy3 = sim_plan_new(method_sp = method_sp, spread = spread, lims = lims, radius = sp_radius,container = "barge")

# Hopper car
sp_xy4 = sim_plan_new(method_sp = method_sp, spread = spread, lims = lims, radius = sp_radius,container = "hopper", compartment = 2, type = "open_top")
sp_xy5 = sim_plan_new(method_sp = method_sp, spread = spread, lims = lims, radius = sp_radius,container = "hopper", compartment = 3, type = "trough")
```

# Results

1. Overlay the sampling plan on the simulated contamination plot.

```{r}
# Discrete spread in 3D space with a 2D projection

# Truck
overlay_draw(method = method_sp, data = rbind(contam_xy, sp_xy1), spread = "discrete", xlim = lims$xlim, ylim = lims$ylim)
overlay_draw(method = method_sp, data = rbind(contam_xy, sp_xy2), spread = "discrete", xlim = lims$xlim, ylim = lims$ylim)

# Barge
overlay_draw(method = method_sp, data = rbind(contam_xy, sp_xy3), spread = "discrete", xlim = lims$xlim, ylim = lims$ylim)

# Hopper car
overlay_draw(method = method_sp, data = rbind(contam_xy, sp_xy4), spread = "discrete", xlim = lims$xlim, ylim = lims$ylim)
overlay_draw(method = method_sp, data = rbind(contam_xy, sp_xy5), spread = "discrete", xlim = lims$xlim, ylim = lims$ylim)
```

```{r, eval = FALSE}
plot_ly(data = subset(contam_sp_xy, subset = label != "sample point"), 
            x = ~X, y = ~Y, z = ~Z, name = "contamination", type = "scatter3d", mode = "markers") %>%
  add_trace(data = subset(contam_sp_xy, subset = label == "sample point"), x = ~X, y = ~Y, z = ~Z, name = "sample", mode = "markers")
```

```{r, eval = FALSE, echo = FALSE}
knitr::include_graphics(path = "example1.png")
```

```{r, echo= FALSE, eval = FALSE}
contam_level_draw(dimension = "3d", method = fun, spread_radius = spread_radius, LOC = LOC, df_contam = contam_xy, xlim = x_lim, ylim = y_lim, interactive = FALSE)
```