annulus_demo

library(ashapesampler)
library(alphahull)
library(ggplot2)
library(doParallel)
#> Loading required package: foreach
#> Loading required package: iterators
#> Loading required package: parallel
library(parallel)
cores <- min(2L, detectCores())

In this document, we demonstrate the \(\alpha\)-shape sampler pipeline by simulating the process of learning a set of two-dimensional shapes (in this case, annuli) and simulating a new shape from that. This vignette requires the packages alphahull, ggplot2, parallel, and doParallel in addition to ashapesampler.

We begin by setting the parameters for our simulation. We will fix \(\alpha=0.25\) and \(n=100\), and draw 20 shapes for our data set. Our true underlying manifold will be the annulus with outer radius 0.75 and inner radius 0.25.

set.seed(201723)
my_alpha = 0.15
n = 1000
N= 10
r_maj <- 0.75
r_min <- 0.25

Next we will draw the shapes themselves.

ann_list <- list()
complex_list <- list()
tau_vec <- vector("numeric", N)
for (k in 1:N){
  ann_pts <- runif_annulus(n, r_maj, r_min)
  ann_list[[k]] <- ashape(ann_pts, alpha = my_alpha)
  complex_list[[k]] <- get_alpha_complex(ann_pts, my_alpha)
  tau_vec[k] <- tau_bound(ann_list[[k]]$x, complex_list[[k]])
}

Now that we have the shapes generated and imported, we want to sample point clouds to combine. We’ll choose 2.

choose_2 <- sample(N,2)
point_cloud = rbind(ann_list[[choose_2[1]]]$x, ann_list[[choose_2[[2]]]]$x)

Then we will have our \(\tau\) bound be a summary statistic of the \(\tau\) found for each input shape. Here, we will use mean, but one can tweak this to see different results. Note that if \(\tau\) is too small, then the random walk won’t be able to execute around the point cloud, but if \(\tau\) is too big, then we risk losing geometric and topological information in the reconstruction.

tau_vec2 = c(tau_vec[choose_2[1]], tau_vec[choose_2[2]])

Now we can take the parameters and generate a new shape and plot it. Here, we assume k_min=2 as we are in two dimensions.

new_annulus <- generate_ashape2d(point_cloud, J=2, tau=min(tau_vec2), 
                                 cores=cores)
#> [1] "Acceptance Rate is 0.9384"


tri_keep = new_annulus$delvor.obj$tri.obj$trlist[which(new_annulus$delvor.obj$tri.obj$cclist[,3]<new_annulus$alpha), 1:3]
dim_tri = dim(tri_keep)[1]
tri_keep = as.vector(t(tri_keep))
triangles = data.frame("id"=sort(rep(1:dim_tri, 3)), 
                       "x"=new_annulus$x[tri_keep, 1], 
                       "y"=new_annulus$x[tri_keep,2])
extremes = as.data.frame(new_annulus$x[new_annulus$alpha.extremes,])

edges = as.data.frame(new_annulus$edges[,3:6])

ggplot(data.frame(new_annulus$x), aes(x=X1, y=X2)) +
  geom_polygon(data=triangles, aes(x=x, y=y, group=id), fill="gray") +
  geom_segment(data=edges, aes(x=x1, y=y1, xend=x2, yend=y2), color="blue")+
  geom_point(data=extremes, aes(x=V1, y=V2), size=1.5)+
  theme_classic()