BARS Transmission Model: Chicago

– Aditya, Nick, Jonathan, and Nikki

Today, we will see results from the Chicago model, for a 25,000 day simulation (approx. 68 years).

The model includes:

• Demography: Arrivals, departures, aging;
• Biology: Individual-level CD4 and viral load trajectories, with and without ART;
• Treatment: Current implementation of ART initiation is simplistic, we are working on more realisic model (details follow);
• Sexual behavior and networks: Including main and casual partnerships, and concurrency, duration and density for each partnership type.

Currently, we are modeling a testing-and-diagnosis piece. The implementation for the starting network is complete. Temporal simulation of testing and diagnosis is in progress. More details are here.

Following this step, we will include PrEP.

Other points for discussion: Lag time between diagnosis and ART initiation; ART/PrEP cessation; mixing parameters; matching simulated and empirical age structure.

Below, is an “interrogation” of the model, to demonstrate its various features.

These datasets include the sexual network at the last time step, biomarker data with detailed trajectories, infection, death and partnership events, and counts of various quantities at each time step.

Results

Demography

We first compute the annual growth rate.

The final population size is 7505, corresponding to a growth rate of 1.49%. (The rate of entry of new individuals is a free parameter, and will be adjusted to reflect data form the three cities).

The age distribution at the end of the simulation is below.

age <- net%v%"age"
plot(
  qplot(age, geom="histogram", binwidth=5, 
       col=I("red"), 
       alpha=I(.2))
  )

Prevalence

The overall prevalence at the end of the simulation was about 46.7%. HIV prevalence as a function of time is shown below.

Incidence

The annual incidence rates over the course of the simulation are shown below.

The mean annual incidence, averaged over the length of the simulation, is 4.66%.

Sexual Networks

Main Partnerships

Mean degree.

The target main mean degree is 0.46 and a plot of the simulated mean degree, over the course of the simulation, is below.

Degree distribution.

The simulated and target momentary degree distributions are below.

Partnership duration.

## Loading required package: data.table

## 
## Attaching package: 'data.table'

## The following objects are masked from 'package:dplyr':
## 
##     between, last

The target mean duration is 911. In the simulated data above, the interquartile range is 224, 1075.

Casual Partnerships

Mean degree.

The target casual mean degree is 0.59 and a plot of the simulated mean degree, over the course of the simulation, is below.

Degree distribution.

The simulated and target momentary degree distributions are below.

## degree0 degree1 degree2 degree3 degree4 degree5 degree6 
##    4083    2641     596     155      21       8       1

Partnership duration.

pt_data_casual <- pt_data[network_type == 1]

# group by p1, p2
pt_data_casual[,.(counts = .N),by=.(p1,p2),][,unique(counts)] # 1 2 4 3

## [1] 1 2 4 3

pt_data_casual[,.(counts = .N),by=.(p1,p2),][,.N,by=counts]

##    counts      N
## 1:      1   8830
## 2:      2 121934
## 3:      4    160
## 4:      3     16

# so focus on just 2 is justified
    ## (old DONOTUSE: res <- pt_data_casual[,.(tick, counts = .N),by=.(p1,p2),][counts == 2,abs(diff(tick))])

res <- pt_data_casual[,.(tick, counts = .N),by=.(p1,p2),][counts == 2,.(dur = abs(diff(tick))), by=.(p1,p2)][,dur]

summary(res)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##     1.0    91.0   219.0   314.7   436.0  3582.0

res_dt <- as.data.frame(res)
  qplot(res, geom="histogram", binwidth=100, 
       col=I("red"), 
       alpha=I(.2))

sum.res <- summary(res)

The target mean duration is 335. In the simulated data above, the interquartile range is 91, 436.

Overlap between main and casual partnerships.

Total mean degree: main and casual partnerships combined.

ART metrics

The ART initiation portion was simplified so that all infected's initate ART 1 year after seroconversion. At presetnt, the proportion of the infected on ART was

  on.art <- which(net%v%"art.status" == 1)
  length(on.art)

[1] 2582

  length(on.art)/length(infected)

[1] 0.7364518

0.74. The total proportion of all people in the population who were on ART is 34.4%.

Biomarkers

   ## Selected upto n.to.follow infecteds at random
   n_to_follow <- 10
   infectees <- inf_event_data$infectee
   uniq_biom_pid <- unique(biom_data$p_id)

   infectee_uniq_biom_pid <- uniq_biom_pid[which(uniq_biom_pid %in% infectees)]
   new <- biom_data[which(biom_data$p_id == infectee_uniq_biom_pid[1]),]

   if (length(infectee_uniq_biom_pid > n_to_follow)){
     #infectee_uniq_biom_pid <- infectee_uniq_biom_pid[1:n_to_follow]#first 10
     infectee_uniq_biom_pid <- sample(infectee_uniq_biom_pid, 
                                      n_to_follow, replace=FALSE)# random 10
   }

   for (i in 2:length(infectee_uniq_biom_pid)){
     new_entry <- biom_data[which(biom_data$p_id == infectee_uniq_biom_pid[i]),]
     new <- rbind(new, new_entry)
   }

   new$p_id <- as.character(new$p_id)

   par(mfrow=c(2,1))
   ggplot(new, aes(x=tick/365, y=viral_load, color=p_id))+
          geom_line()+xlab("year")

   ggplot(new, aes(x=tick/365, y=cd4_count, color=p_id))+
          geom_line()+xlab("year")

Example CD4 and viral load trajectories for 10 randomly selected, HIV-infected individuals are shown above. Individuals on ART are identifiable as those whose CD4 counts start to recover after decreasing initially; individuals who do not go on ART do experience a monotonic decrease.

Conclusion

The features we have in the model appear to work well so far. We are in the process of adding more (as mentioned above), and continuing to monitor emergent trends in the behavior of key parameters and outcomes.

Suggestions welcome!