Creating shiny app for spatial disease modelling through time

Hi,
Apologies for the newbie Shiny question....

My basic question is... Does shiny have to have a static dataframe to work from?!

I am trying to produce a shiny app which displays a contact network which changes through time. When I do this in rstudio I just make the contact netwok plot update as the model runs so the nodes appear to move and also change colour depending on their infection status through time. But when I transfer the code to shiny it just waits for the entire model to run and then displays the final time point plot.
I have been experimenting with various ways around this such as creating animated plots in gganimate but again it is very slow to run and I havent managed to fully implement it in shiny yet nut I am wondering if I am missing a really simple fix for this issue?!
Many thanks in advance,
Dave

Hi Dave,
I don't know if I will be able to help, but if you include your code it will be easier for everyone here to help you :slight_smile:

And maybe this can be of help too if you hadn't read it before: https://shiny.rstudio.com/articles/reactivity-overview.html

2 Likes

Hi, and welcome!

Please see the FAQ: What's a reproducible example (`reprex`) and how do I do one? Using a reprex, complete with representative data will attract quicker and more answers.

Dave, having to reverse engineer a problem from a general description can be a big deterrent. Even someone like me, who is not a bright shiny object may be able to help out. But otherwise, the best that can be hoped is to attract the eye of someone who's had the similar challenge, and that requires being lucky with the question title.

1 Like

Thanks, sorry I should include the code I am trying to 'Shiny'! Please see below. The model outputs two different plots that are updated as the model runs which is simulating the passage of time. I am trying to convert this to a shiny app so the user can adjust a couple of the parameters and then see the impact they have on disease spread. I have got some simpler versions of the code to work but I am struggling to get the code to display the plots as they transition through time... it seems that Shiny runs the code and allows the model to finish and then just plots the final situation. I am experimenting with the reactive elements but not quite sure if this is something that can be done... can the plot continually update as the model runs?? I have tried to neaten up the code to a series of functions and then one which 'steps' the process forward one time step but I am struggling to get it to plot one stage at a time and continually update as such. The code I have added runs fine as stand alone... outside of shiny. But if anyone could point me in the right direction as to what functions/commands in Shiny I need to be using I would be very grateful. Many thanks.

########Create and plot associations of Simuland citizens who choose social contacts
########according to degree (number of citizens contacted per day) and loyalty (scales
########the size of the neighborhood explored by each citizen)

rm(list=ls())
set.seed<-1

library(igraph)
par(mfrow=c(2,1))

businessarea=0.1
degree<-5
loyalty<-0.9

#set colours for plotting
cols<-c("white","yellow","orange","red","purple","cyan","black")

pop_set_up<-function(family_size=5,households=100,xMin=-.5,xMax=.5,yMin=-.5,yMax=.5,businessarea){
  
  pop<-family_size*households
  
  #Extended simulation windows parameters
  rExt=businessarea; #extension parameter -- use cluster radius
  xMinExt=xMin-rExt;
  xMaxExt=xMax+rExt;
  yMinExt=yMin-rExt;
  yMaxExt=yMax+rExt;
  #rectangle dimensions
  xDeltaExt=xMaxExt-xMinExt;
  yDeltaExt=yMaxExt-yMinExt;
  areaTotalExt=xDeltaExt*yDeltaExt; #area of extended rectangle
  
  #homestead positions
  #x and y (uniform) coordinates of Poisson points for the parent
  xParent=xMinExt+xDeltaExt*runif(households);
  yParent=yMinExt+yDeltaExt*runif(households);
  
  #replicate for daughters 
  #replicate parent points (ie centres of disks/clusters)
  xP=rep(xParent,each=family_size);
  yP=rep(yParent,each=family_size);
  
  ###from epidemic model...setup
  #age the citizens
  age<-rbinom(length(xP),1,0.25)
  
  output<-list(pop,xP,yP,age,households,xMin,xMax,yMin,yMax)
  names(output)<-c("pop","xP","yP","age","households","xMin","xMax","yMin","yMax")
  return(output)
  
}

disease_set_up<-function(pop){
  
  #start with 2 infected individuals
  expa<-sample(1:pop,2,replace=FALSE)
  
  #a necessary disease parameter
  E_I1<-5
  
  #create dataframe to store disease state
  S<-rep(1,pop)
  E<-rep(0,pop)
  I1<-rep(0,pop)
  I2<-rep(0,pop)
  I3<-rep(0,pop)
  R<-rep(0,pop)
  D<-rep(0,pop)
  status<-data.frame(S,E,I1,I2,I3,R,D)
  
  status$S[expa]<-0
  status$E[expa]<-1
  
  d_exp<-rep(NA,pop)
  d_inf1<-rep(NA,pop)
  d_inf2<-rep(NA,pop)
  d_inf3<-rep(NA,pop)
  
  d_exp[status$E==1]<-rpois(1,E_I1)
  
  #create matrix to store epidemic progress
  progression<-matrix(0,nr=1,nc=ncol(status))
  progression[1,]<-colSums(status)
  
  output<-list(status,progression,d_exp,d_inf1,d_inf2,d_inf3)
  names(output)<-c("status","progression","d_exp","d_inf1","d_inf2","d_inf3")
  return(output)
  
}


disease_model<-function(Spop=pop$pop,mat,age=pop$age,status=dis$status,d_exp=dis$d_exp,d_inf1=dis$d_inf1,d_inf2=dis$d_inf2,d_inf3=dis$d_inf3,progression=dis$progression){
  #Probability of becoming exposed having contacted an infectious individual (daily)
  #Need to work out R0 based on other parameters
  S_E<-0.1
  #lambda for a Poisson draw for the length of this period
  E_I1<-5
  #probability of transitioning to serious case for young (daily)
  yI1_I2<-0.01
  #and same for old
  oI1_I2<-0.1
  #diff between young and old
  diffI1_I2<-oI1_I2-yI1_I2
  #probability of transitioning to critical case for young (daily)
  yI2_I3<-0.05
  #and same for old
  oI2_I3<-0.2
  #diff between young and old
  diffI2_I3<-oI2_I3-yI2_I3
  #probability of death for young (daily)
  yI3_D<-0.1
  #and same for old
  oI3_D<-0.4
  #diff between young and old
  diffI3_D<-oI3_D-yI3_D
  #lambda for a Poisson draw for duration of a mild infection
  yI1_R<-7
  #and same for old
  oI1_R<-7
  #diff between young and old
  diffI1_R<-oI1_R-yI1_R
  #lambda for a Poisson draw for duration of a hospitalised infection
  yI2_R<-5
  #and same for old
  oI2_R<-5
  #diff between young and old
  diffI2_R<-oI2_R-yI2_R
  #lambda for a Poisson draw for duration of a critical infection
  yI3_R<-5
  #and same for old
  oI3_R<-5
  #diff between young and old
  diffI3_R<-oI3_R-yI3_R
  
  #this will block individuals changing twice in the same time-step
  f2c<-rep(1,Spop)
  inf<-sign(status$I1+status$I2+status$I3)
  risk<-rep(0,Spop)
  risk[as.numeric(names(table(which(inf*mat==1,arr.ind=TRUE)[,2])))]<-table(which(inf*mat==1,arr.ind=TRUE)[,2])
  prob<-(1-S_E)^risk
  infect<-rbinom(Spop,1,1-prob)*status$S
  status$S[which(infect==1)]<-0
  status$E[which(infect==1)]<-1
  d_exp[!is.na(d_exp)&d_exp>0]<-d_exp[!is.na(d_exp)&d_exp>0]-1
  d_exp[which(infect==1)]<-rpois(length(which(infect==1)),E_I1)
  f2c[which(infect==1)]<-0
  EI1<-which(status$E==1&!is.na(d_exp)&d_exp==0)
  status$E[EI1]<-0
  status$I1[EI1]<-1
  f2c[EI1]<-0
  I1I2<-rbinom(Spop,1,yI1_I2+age*diffI1_I2)*status$I1*f2c
  status$I1[which(I1I2==1)]<-0
  status$I2[which(I1I2==1)]<-1
  f2c[which(I1I2==1)]<-0
  I2I3<-rbinom(Spop,1,yI2_I3+age*diffI2_I3)*status$I2*f2c
  status$I2[which(I2I3==1)]<-0
  status$I3[which(I2I3==1)]<-1  
  f2c[which(I2I3==1)]<-0
  I3D<-rbinom(Spop,1,yI3_D+age*diffI3_D)*status$I3*f2c
  just.died<-which(I3D==1)
  status$I3[just.died]<-0
  status$D[just.died]<-1
  f2c[just.died]<-0
  
  EI1y<-EI1[EI1%in%which(age==0)]
  EI1o<-EI1[EI1%in%which(age==1)]
  d_inf1[!is.na(d_inf1)&d_inf1>0]<-d_inf1[!is.na(d_inf1)&d_inf1>0]-1
  if(length(EI1y)>0){
    d_inf1[EI1y]<-rpois(length(EI1y),yI1_R)
  }
  if(length(EI1o)>0){
    d_inf1[EI1o]<-rpois(length(EI1o),oI1_R)
  }
  
  d_inf2[!is.na(d_inf2)&d_inf2>0]<-d_inf2[!is.na(d_inf2)&d_inf2>0]-1
  d_inf2[which(I1I2==1&age==0)]<-rpois(length(which(I1I2==1&age==0)),yI2_R)
  d_inf2[which(I1I2==1&age==1)]<-rpois(length(which(I1I2==1&age==1)),oI2_R)
  
  d_inf3[!is.na(d_inf3)&d_inf3>0]<-d_inf3[!is.na(d_inf3)&d_inf3>0]-1
  d_inf3[which(I2I3==1&age==0)]<-rpois(length(which(I2I3==1&age==0)),yI3_R)
  d_inf3[which(I2I3==1&age==1)]<-rpois(length(which(I2I3==1&age==1)),oI3_R)
  
  I1R<-which(status$I1==1&!is.na(d_inf1)&d_inf1==0)
  status$I1[I1R]<-0
  status$R[I1R]<-1
  f2c[I1R]<-0
  I2R<-which(status$I2==1&!is.na(d_inf2)&d_inf2==0)
  status$I2[I2R]<-0
  status$R[I2R]<-1
  f2c[I2R]<-0
  I3R<-which(status$I3==1&!is.na(d_inf3)&d_inf3==0)
  status$I3[I3R]<-0
  status$R[I3R]<-1
  f2c[I3R]<-0
  
  progression<-rbind(progression,colSums(status))
  
  status.now<-as.matrix(status)%*%as.matrix(seq(1,7))
  
  output<-list(status,status.now,progression,d_exp,d_inf1,d_inf2,d_inf3,just.died)
  names(output)<-c("status","status.now","progression","d_exp","d_inf1","d_inf2","d_inf3","just.died")
  return(output)
  
}

night_net<-function(pop=pop,households=households){
  
  ##########Create association matrix (Night)...will stay constant unless Sick or dead go to hospital
  #set up association matrix
  ass.night<-array(0,dim=c(pop,pop))
  #a useful index vector
  assseq<-seq(1,pop)
  
  #make association matrix
  #I feel like this could be speeded up, avoiding a for-loop
  #just want bliocks of 1s for housemates
  for(i in 1:households){
  ass.night[(i*5-4):(i*5),(i*5-4):(i*5)]<-1}
  
  return(ass.night)
  
}

forward<-function(Spop=pop$pop,businessarea,xP=pop$xP,yP=pop$yP,loyalty,degree,ass.night,status){
  
  #####DAY
  #Generate the (relative) locations in polar coordinates by 
  #simulating independent variables.
  theta=2*pi*runif(Spop); #angular coordinates 
  rho=businessarea*sqrt(runif(Spop)); #radial coordinates 
  
  #Convert from polar to Cartesian coordinates
  x0=rho*cos(theta);
  y0=rho*sin(theta);
  
  #translate points (ie parents points are the centres of cluster disks)
  x=xP+x0;
  y=yP+y0;
  
  ##########Bit of data handling
  citizens<-data.frame(x,y)
  
  ##########Create association matrix (Business)
  #calculate distance matrix
  dm<-as.matrix(dist(citizens,method="euclidean"))
  
  #set up (empty) association matrix
  ass.day<-dm*0
  
  ###here we loop through each citizen, have them select "degree" partners, and associate them symmetrically.
  ###this is not the strict definition of "degree", but it does seem to describe human behaviour ("I seek n contacts per day")
  #loop through citizens
  for(assrow in 1:Spop){
    dist_t<-dm[assrow,]
    possible.partners<-seq(1,Spop)
    possible.partners<-possible.partners[dist_t<=quantile(dist_t,1-loyalty)&dist_t>0]
    partners<-sample(possible.partners,size=degree,replace=F)
    ass.day[assrow,partners]<-ass.day[assrow,partners]+1
  }
  
  #ass.day<-sapply(1:nrow(ass.day),assign_contacts)
  ass.day<-sign(ass.day+t(ass.day))
  
  ##MJS ADDITION##
  #Disconnect individuals that have severe infection from day time network
  ass.day[status$I2==1,]<-0
  ass.day[,status$I2==1]<-0
  ass.day[status$I3==1,]<-0
  ass.day[,status$I3==1]<-0
  
  #Dave query....doi we want to do something interesting with them like leave them at home so that folk could still interact with them?
  #or have them iunteract with health professionals, hence affecting the health service carrying capacity?
  
  #this bit makes an association network for plotting,removing any dead individs
  business.ass<-ass.day
  business.ass[status$D==1,]<-0
  business.ass[,status$D==1]<-0
  
  mat<-ass.night+ass.day
  
  #simplify the total association matrix
  mat<-sign(mat)
  #remove any deaths from the total association matrix
  mat[status$D==1,]<-0
  mat[,status$D==1]<-0

  dis<-disease_model(mat=mat)
  
  output<-list(business.ass,mat,dis,citizens)
  names(output)<-c("business.ass","mat","dis","citizens")
  return(output)
    
}

####################################
####################################
####################################

pop<-pop_set_up(businessarea=businessarea)
dis<-disease_set_up(pop=pop$pop)
an<-night_net(pop=pop$pop,households=pop$households)
status<-dis$status
counter<-0
t<-0
while((counter<10)&(t<200)){
  t<-t+1
  step<-forward(businessarea=businessarea,degree=degree,loyalty=loyalty,ass.night=an,status=status)
  dis<-step$dis
  status<-dis$status
  if(sum(status$R)+sum(status$D)==nrow(status)){counter<-counter+1}
  print(dis$progression[nrow(dis$progression),])

#plot and add association lines
#cool code finds the array indices of the interacting dyads
dd<-which(step$business.ass>0,arr.ind=T)
#set up array of arrow vertices
arrow.x1<-arrow.x2<-arrow.y1<-arrow.y2<-numeric(dim(dd)[1])
#and for each dyad, populate the vertices. Might be able to avoid a for-loop here?
for(i in 1:dim(dd)[1]){
arrow.x1[i]<-step$citizens$x[dd[i,1]]
arrow.y1[i]<-step$citizens$y[dd[i,1]]
arrow.x2[i]<-step$citizens$x[dd[i,2]]
arrow.y2[i]<-step$citizens$y[dd[i,2]]}
#and for plotting purposes, checks whether alive
is.alive<-dis$status$D==0
is.ill<-(dis$status$I2==1|dis$status$I3==1)
step$x[is.ill]<-pop$xP[is.ill]
step$y[is.ill]<-pop$yP[is.ill]

dev.hold()
plot(NULL,xlim=c(pop$xMin,pop$xMax),ylim=c(pop$yMin,pop$yMax),xlab="latitude",ylab="longitude",pch=21,bg=cols[(dis$status.now[is.alive])],col="black",cex=1.5)
arrows(arrow.x1,arrow.y1,arrow.x2,arrow.y2,code=3,length=0,col="gray80")
points(step$citizens$x[is.alive],step$citizens$y[is.alive],pch=21,bg=cols[(dis$status.now[is.alive])],col="black",cex=1.5)
points(pop$xP[dis$just.died],pop$yP[dis$just.died],pch=21,bg="black",col="black",cex=1.5)

severe.cases<-(dis$progression[,4]+dis$progression[,5])/pop$pop
mort<-dis$progression[,7]/pop$pop
Time<-seq(1,201)
plot(NULL,xlab="Time",ylab="Cumulative mortality",ylim=c(0,0.1),xlim=c(1,Time[t+1]),type="n",las=1)
lines(mort[1:(t+1)]~Time[1:(t+1)],lwd=2)
lines(severe.cases[1:(t+1)]~Time[1:(t+1)],ylab="Severity",lwd=2,col="red")
dev.flush()


}
  
1 Like

Take a look at the gganimate package, which is designed to step through a series of evolving plot objects. You may want to consider dividing the plots, though, since an animated 2x2 window may be difficult to assimilate.

Thanks technorat.... yeah I have been trying to play with the gganimate package... I have made some progress on the coding today so I am optimistic!!
Cheers,
Dave

1 Like