Local curvature

Here we introduce the definition of the local curvature. In two-dimensional space, local curvature is defined by \[c = \frac{|x^{\prime}y^{\prime\prime} - x^{\prime\prime}y^{\prime}|}{\sqrt{[(x^{\prime})^2 + (y^{\prime})^2]^3}},\] where \(x^{\prime}\) and \(x^{\prime\prime}\) are the first and second derivatives, respectively. In shape analysis, the coordinates are rotated so that the downward direction in Eucledean space is toward the center of the debris when calculating the local curvature at every pixel. Local curvature values are obtained at each pixel around the boundary. Negative curvature value at a boundary pixel indicates that the boundary is bent outward (e.g., circle has negative curvature throughout it) at the pixel, whereas positive one means that it is bent inward (i.e., dent). Hence frequent changes in curvature over the boundary can be interpreted to be more “wiggly”.

Hypotheses test

It is difficult to make a typical assumptions commonly used in statistical hypotheses such normality or indepence sampling in our setting, because of the spatial dependence inherent in microscopic images. Thus we consider a nonparametric test that requires minimal assumptions and conditions. Our test relies on resampling procedure. For simplicity, we assume that same number of sample images are obtained for both plastic types, although it can be flexibly chosen with minimal modification. The procedure is summarized below.

  1. For each sample (\(i=1, \ldots, I\)) and type of the plastic \(j\) = 1(ABS), 2(PVC), randomly pick a partition of length \(K\). Denote the partition \(y_{ij}^{(b)} = ( y_{ij1}^{(b)}, y_{ij2}^{(b)}, \ldots, y_{ijK}^{(b)} )\).
  2. Within this partition, obtain the length \(K-1\)-vector containing the increment (or decrement) at each point, denoted by \(\tilde{y}_{i1}^{(b)}\) and \(\tilde{y}_{i2}^{(b)}\), respectively.
  3. Obtain the average curvature change for \(j\), denoted by \(\bar{y}_{j}^{(b)} = \frac{1}{I} \frac{1}{K-1} \sum_{i=1}^{I} \sum_{k=1}^{K-1} |\tilde{y}_{ijk}^{(b)}|\), and define \(\bar{w}^{(b)} = (\bar{y}_{1}^{(b)}, \bar{y}_{2}^{(b)})\).

After repeating 1-3 for \(b = 1,\ldots, B\), \(\bar{w}^{(1)}, \ldots, \bar{w}^{(B)}\) is available. These statistics can be used for other inferences, such as conducting a hypothesis test or constructing a confidence interval. In this paper, we focus on the comparison of the average curvature for different plastic types. Specifically, we conduct a hypothesis test with the following hypotheses:

Apparently this test is not available from the “over-the-shelf” statistical software. We made an R package named microplastics for analysis. The package can be installed by running

install.packages("devtools")
devtools::install_github("ydhwang/microplastics")

and can be used by

library(microplastics)

Installation is needed only once when used for the first time.

Image loading, preparation, and hypothesis test

We assume that the image files have png extension and are stored in separate directories named ../piece/pvc and ../piece/abs according to their types. They are stored into two list named PVC_list and ABS_list.

PVC_list <- list.files(path="../piece/pvc", pattern = "*.png") %>%  str_to_lower()
PVC <- list()
PVC_nms <- PVC_list %>% str_remove(".png")
for (k in 1:length(PVC_list)){
  PVC[[PVC_nms[k]]] <- OpenImageR::readImage(paste0("../piece/pvc/", PVC_list[k]))
}

ABS_list <- list.files(path="../piece/abs", pattern = "*.png") %>%  str_to_lower()
ABS <- list()
ABS_nms <- ABS_list %>% str_remove(".png")
for (k in 1:length(ABS_list)){
  ABS[[ABS_nms[k]]] <- OpenImageR::readImage(paste0("../piece/abs/", ABS_list[k]))
}

This is the original image:

knitr::include_graphics("../piece/abs/abs_1.png")

Once the images are loaded as a matrix, they can be seen (with some rotation):

image(ABS[[1]][,,3], col = terrain.colors(100), axes = FALSE)

Now we proceed with some steps to adjust the image size and covert the image to matrices.

ABS_m <- list()
for (k in 1:4){ 
  ABS_m[[k]] <- EBImage::resize(ABS[[k]][,,4], 200, 200)
}

PVC_m <- list()
for (k in 1:4){ 
  PVC_m[[k]] <- EBImage::resize(PVC[[k]][,,4], 200, 200)
}

Now the estimated curvature along the boundaries is stored in abs_cur_val and pvc_cur_val.

abs_cur_val <- list()
for (k in seq_along(ABS_m)){
  abs_cur_val[[k]] <- get_curvature(ABS_m[[k]], 10)
}

pvc_cur_val <- list()
for (k in seq_along(PVC_m)){
 pvc_cur_val[[k]] <- get_curvature(PVC_m[[k]], 10)
}

This is an example of local curvature.

get_curvature_image(ABS_m[[1]], 10)

Now we flatten each sample into a series of local curvature estimates. Local curvatures are estimated for each plastic type, which are depicted in the figure below. id represents the index for the boundary pixels.

(curv_figure <- 
rbind(df_by_sample(abs_cur_val, 'ABS'),
      df_by_sample(pvc_cur_val, 'PVC')) %>% 
  ggplot() + aes(x = id, y = curvature) + 
  geom_line(lwd = 1) + scale_x_continuous() + theme_bw() +
  facet_grid(sample_id ~ key) + ylab("Curvature") + xlab("ID"))
Local curvature estimates for plastic types and samples

Local curvature estimates for plastic types and samples

One can see the PVC has more abrupt change points, which indicates that it has more wiggly boundaries.

Now we reshape the images into tibble structure that can be handy for statistical shape analysis.

sample_abs <- df_by_sample(abs_cur_val, 'ABS')
sample_pvc <- df_by_sample(pvc_cur_val, 'PVC')
sample_combined <- rbind(sample_abs, sample_pvc)

Set the parameters for sampling and conduct a hypothesis test.

K <- 100 # length K partition
B <- 1000

# B = 1000 bootstrap test
boot <- list()

for (b in 1:B){
  boot[[b]] <- perm_test(sample_combined)
  if (b%%100 == 0) cat(b, "\n")
}
## 100 
## 200 
## 300 
## 400 
## 500 
## 600 
## 700 
## 800 
## 900 
## 1000
boot_tbl <- bind_rows(boot) %>% rename(Index = key)

boot_tbl %>% 
  ggplot() + aes(fill = Index, x = M) + geom_density(alpha = 0.2) + xlab("Curvature Change") + ylab("Density")