Skip to main content

Markov modelling of viral load adjusting for CD4 orthogonal variable and multivariate conditional autoregressive mapping of the HIV immunological outcomes among ART patients in Zimbabwe

Abstract

Background

This study aimed to jointly model HIV disease progression patterns based on viral load (VL) among adult ART patients adjusting for the time-varying “incremental transients states” variable, and the CD4 cell counts orthogonal variable in a single 5-stage time-homogenous multistate Markov model. We further jointly mapped the relative risks of HIV disease progression outcomes (detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) at the last observed visit) conditional not to have died or become loss to follow-up (LTFU).

Methods

Secondary data analysis of individual-level patients on ART was performed. Adjusted transition intensities, hazard ratios (HR) and regression coefficients were estimated from the joint multistate model of VL and CD4 cell counts. The mortality and LTFU transition rates defined the extent of patients’ retention in care. Joint mapping of HIV disease progression outcomes after ART initiation was done using the Bayesian intrinsic Multivariate Conditional Autoregressive prior model.

Results

The viral rebound from the undetectable state was 1.78times more likely compared to viral suppression among patients with VL ranging from 50-1000copies/uL. Patients with CD4 cell counts lower than expected had a higher risk of viral increase above 1000copies/uL and death if their VL was above 1000copies/uL (state 2 to 3 (λ23): HR = 1.83 and (λ34): HR = 1.42 respectively). Regarding the time-varying effects of CD4 cell counts on the VL transition rates, as the VL increased, (λ12 and λ23) the transition rates increased with a decrease in the CD4 cell counts over time. Regardless of the individual’s VL, the transition rates to become LTFU decreased with a decrease in CD4 cell counts. We observed a strong shared geographical pattern of 66% spatial correlation between the relative risks of detectable VL and immune deterioration after ART initiation, mainly in Matabeleland North.

Conclusion

With high rates of viral rebound, interventions which encourage ART adherence and continual educational support on the barriers to ART uptake are crucial to achieve and sustain viral suppression to undetectable levels. Area-specific interventions which focus on early ART screening through self-testing, behavioural change campaigns and social support strategies should be strengthened in heavily burdened regions to sustain the undetectable VL. Sustaining undetectable VL lowers HIV transmission in the general population and this is a step towards achieving zero HIV incidences by 2030.

Background

Antiretroviral therapy (ART) has since been the backbone of HIV prevention and control as it reduces viral load (VL) replication in the human host by blocking the virus life cycle [1]. Once the virus replication has been inhibited, the CD4 cell counts increase, and the individual life expectancy is expanded [2]. HIV patients’ management involves monitoring VL and CD4 cell counts prognostic markers through laboratory repeated measurements. The immunological markers can be utilized in understanding the HIV disease progression patterns among ART patients [2].

Multistate Markov models are mathematical models which have been used to evaluate HIV disease processes; however, the number of states, state cuts-off points and the number of transitions vary across studies [3, 4]. This explains the flexibility of the multistate models’ implementation, but the models become complex as the number of states and transitions increase. These models have been used extensively in HIV disease progression using CD4 cell counts states [5, 6] or VL states [7, 8] to define the model states (categories) separately. Of these two prognostic markers, VL is the preferred marker in HIV monitoring due to its high sensitivity. However, there had been a delay in the implementation and rolling out of VL testing in most developing countries due to cost-related challenges [9].

The low-middle-income countries (LMIC) have traditionally relied on the use of CD4 cell counts in HIV disease progression monitoring as this has been the readily available laboratory marker; however, if both VL and CD4 cell counts are available, it is imperative to model these prognostic markers jointly to understand better the HIV disease progression patterns. Joint modelling of these two prognostic markers helps to explain those effects which one marker cannot explain in the absence of the other marker.

In HIV disease monitoring programmes, Bayesian spatial modelling is an emerging tool to analyze spatially related multidimensional data with an underlying spatial process to guide policy [8, 10]; however, joint spatial modelling using the Bayesian intrinsic Multivariate Conditional Autoregressive (MCAR) prior has not been fully utilized in this field. The advantage of joint mapping is that it gives an understanding of the HIV dynamics and the spatial overlap between the joint mapped outcomes.

The main objective of our study was to jointly model HIV disease progression using two 5-stage time-homogenous multistate Markov models based on CD4 cell counts and VL. We firstly fitted two time-homogeneous multistate Markov models with states defined by CD4 cell counts and VL. In each multistate model, we jointly model these prognostic markers with one marker defining the finite multistate model states and the other marker forming the covariate matrix components of the regression model. The first covariate was an orthogonal variable generated using the principal component analysis (PCA) [11]. The second covariate was the time-varying “incremental transient states” variable to estimate the changes in transition rates over time with respect to how the marker changes [7]. Uniquely to this study is the inclusion of these two covariates in a single multistate model covariate matrix which has been a gap in previous studies that incorporated either the orthogonal variable or the time-varying covariate only [11]. We further performed a joint spatial mapping of HIV disease progression outcomes (detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) at the last observed visit conditional not to have died or become loss to follow-up (LTFU) to describe the spatial overlap of the two outcomes.

Methods

Data source, description and study design

This study was a secondary analysis of data from patients’ records compiled for monitoring and guiding programme planning. The data used in this study came from the Zimbabwe National ART programme collected through the electronic patients’ management system (ePMS) database described elsewhere [12, 13]. The electronic system was implemented nationally to improve and increase efficiency in HIV patients’ management and monitor their response to ART. A stratified sampling of health facilities offering HIV services was done; hence, a representative sample of primary, secondary, tertiary and quaternary health facilities which provide HIV services was achieved during the ePMS roll out. However, a continual up-scale is currently ongoing, which will ensure that the system covers all health facilities.

We considered individuals aged 15 years and above, with at least two repeated measurements of both CD4 cell counts and VL who initiated ART between 2004 and 2017. However, since this was programme data, the follow-up visits of the patients were intermittent, and there may have been clinical reasons for requesting CD4 count and VL measurements from patients. Each individual had an average follow-up period of 2 years, and the majority of the patients (79.2%) in the final sample considered for analysis were enrolled on ART between 2009 and 2015.

The time-homogeneous Markov multistate model formulation

A multistate model is a stochastic continuous-time process {X(t), t [0, T)years} defined as a finite space X = {1, 2, 3, 4, 5} based on VL or CD4 cell counts states:

$$ \mathrm{Viral}\ \mathrm{load}\ \left(\mathrm{VL}\right)\ \mathrm{states}\ \left(\mathrm{copies}/\mathrm{uL}\right)=\left\{\begin{array}{l}1;\mathrm{VL}<50\\ {}2;50\le \mathrm{VL}<1000\\ {}3;\mathrm{VL}\ge 1000\\ {}4;\mathrm{Dead}\\ {}5;\mathrm{Loss}\ \mathrm{to}\ \mathrm{follow}\hbox{-} \mathrm{up}\end{array}\right. $$
(1)
$$ \mathrm{CD}4\ \mathrm{cell}\ \mathrm{counts}\ \mathrm{states}\ \left(\mathrm{cells}/\mathrm{uL}\right)=\left\{\begin{array}{l}1;\mathrm{CD}4\ge 500\\ {}2;350\le \mathrm{CD}4<500\\ {}3;\mathrm{CD}4<350\ \\ {}4;\mathrm{Dead}\\ {}5;\mathrm{Loss}\ \mathrm{to}\ \mathrm{follow}\hbox{-} \mathrm{up}\end{array}\right. $$
(2)

where LTFU was defined as a failure of a patient to report for drug refill for at least 90 days from the last appointment date or if the patient missed the next scheduled visit date and never showed up again. The schematic presentation of the 5-stage multistate model is shown in Fig. 1. Based on these possible transitions shown in Fig. 1, the corresponding transition rates are defined by a 5 × 5 transition matrix Q(t) with λij elements defining the movement between state i and state j with properties \( {\sum}_{j=1}^5{\lambda}_{ij}=0 \) and λii =  − ∑i ≠ jλij. In reality, the individuals who become LTFU can return to the clinic for continual monitoring. However, in this study, we could not ascertain any returns of these participants after they became LTFU. This means the λ5j transition rates from the LTFU state were not estimated as the data could not support this.

$$ Q(t)=\left(\begin{array}{ccccc}-\left({\lambda}_{12}+{\lambda}_{13}+{\lambda}_{14}+{\lambda}_{15}\right)& {\lambda}_{12}& {\lambda}_{13}& {\lambda}_{14}& {\lambda}_{15}\\ {}{\lambda}_{21}& -\left({\lambda}_{21}+{\lambda}_{23}+{\lambda}_{24}+{\lambda}_{25}\right)& {\lambda}_{23}& {\lambda}_{24}& {\lambda}_{25}\\ {}{\lambda}_{31}& {\lambda}_{32}& -\left({\lambda}_{31}+{\lambda}_{32}+{\lambda}_{34}+{\lambda}_{35}\right)& {\lambda}_{34}& {\lambda}_{35}\\ {}0& 0& 0& 0& 0\\ {}{\lambda}_{51}& {\lambda}_{52}& {\lambda}_{53}& {\lambda}_{54}& -\left({\lambda}_{51}+{\lambda}_{52}+{\lambda}_{53}+{\lambda}_{54}\right)\end{array}\right) $$
Fig. 1
figure 1

The schematic diagram for the possible transition between the defined state for both the viral load and the CD4 cell counts models [State 1(VL < 50or CD4 ≥ 500), State 2(50 ≤ VL < 1000 or 350 ≤ CD4 < 500), State 3(VL ≥ 1000 or CD4 < 350), State 4 (dead) and State 5 (Loss to follow-up (LTFU)]

The effects of the covariates on each model were modelled using the semi-parametric proportionality hazards. We aimed to describe how the uncorrelated factor of VL values can explain the component of mortality, LTFU or HIV disease progression transition rates that cannot be explained by the CD4 cell counts alone ignoring VL measurements. To achieve this, we generated an orthogonal variable through the PCA technique, which is a data reduction approach used to combine highly correlated variables into uncorrelated components to improve model efficiency as described elsewhere [11]. The orthogonal variable was then included as a covariate in the proportionality hazard model. The second covariate was the time-varying VL measurements variable which was categorized into “incremental transient states”. The cut-off points for the time-varying VL covariate were similar to those defined in equation [2] above, excluding the mortality and the LTFU states. We simultaneously modelled the effects of these two covariates in a single model.

Therefore, we fitted two proportional hazard multistate models of the generic form:

$$ {\lambda}_{ij,k}/\boldsymbol{Z}(t)={\lambda}_{ij,(0)}\exp \left({\beta}_{ij}^{\prime }{\boldsymbol{Z}}_k\right),i\ne j $$
(3)

where λij, k/Z(t) is the transition rate between state i and state j at the time t given a covariate matrix Z and λij, (0) is the baseline hazard rate of the model. The two multistate Markov models fitted were:

$$ {\lambda}_{ij,k(CD4)}/\boldsymbol{Z}(t)={\lambda}_{ij\left( CD4,0\right)}\exp \left({\beta}_{ij, CD4(1)}{P}_{VL(k)}^{\ast }+{\beta}_{ij, CD4(2)}{VL}_{states,k}\right) $$
(4)
$$ {\lambda}_{ij,k(VL)}/\boldsymbol{Z}(t)={\lambda}_{ij\left( VL,0\right)}\exp \left({\beta}_{ij, VL(1)}{P}_{CD4(k)}^{\ast }+{\beta}_{ij, VL(2)} CD{4}_{states,k}\right) $$
(5)

where equation [4] defines the CD4 cell counts multistate Markov model with λij(CD4, 0) as the baseline transition intensity for an individual k with positive residual values for the orthogonal effect of VL \( \left({P}_{VL(k)}^{\ast }=0\right) \) and the time-varying VL transient levels (VLstates, k = 1(VL < 50 copies/uL)). Similarly, equation [5] defines the VL multistate Markov model with λij(VL, 0) representing the baseline transition intensities for an individual k with positive residual values for the orthogonal effect CD4 cell counts \( \left({P}_{CD4(k)}^{\ast }=0\right) \) and the time-varying CD4 cell counts levels (CD4states, k = 1(CD4 ≥ 500 cells/uL)). The βij, (CD4)d and βij, (VL)d for d = 1, 2 on both the CD4 cell count model and the VL model are the log-linear effects (coefficients) of the corresponding orthogonal and time-varying covariates.

Since the analyses were purely based on existing modelling approached within the “msm” R library, additional covariates like baseline CD4 cell counts, baseline VL values, age, sex and WHO staging could be adjusted for in the model; however, in this study, we encountered a convergence warning as more covariates were added. When the multistate Markov model fails to converge, it means the optimization criterion could not converge to the maximum likelihood; hence, the standard errors of the estimated parameters were not calculated, and the confidence interval estimates were missing on the printout of results. This modelling challenge is normally faced if there is data sparsity in some cells within the multistate model. Usually, the scaling factor is increased to a higher value adjusted to normalize the likelihood and prevent overflow within the optimization process; however, in our case, this was not helpful.

The final model selection was based on the Akaike Information Criterion (AIC) values which are defined as AIC =  − 2(Log − likelihood) + 2q, where 2q represents the variance component, q is the number of parameters to be estimated in the fitted model and the bias is defined by −2(Log − likelihood). The better fitting model is the model with the lowest AIC value. The nested models were assessed using the likelihood ratio test (LRT) defined as:

$$ LRT=-2\times {\log}_e\left({L}_1\left(\hat{\phi}\right)/{L}_2\left(\hat{\phi}\right)\right) $$
(6)

where \( {L}_1\left(\hat{\phi}\right) \) is the likelihood for the simple (unsaturated) model with few covariates and \( {L}_2\left(\hat{\phi}\right) \) is the likelihood for the full (saturated) model with additional covariate(s). A significant p-value < 0.05 for the LRT leads to the rejection of the null hypothesis (simple model is better) in favour of the alternative hypothesis (full model is better). To assess how well the final selected models predict HIV disease, mortality and LTFU, percentage prevalence in each state were plotted to compare the observed and the expected frequencies.

Joint mapping of HIV progression based on two prognostic markers

Furthermore, an intrinsic MCAR prior model by Besag et al. (1991) [14] was fitted to jointly map the relative risks of the two HIV disease progression immunological outcomes (immune deterioration (CD4 < 350cells/uL) and detectable VL (VL≥50copies/uL) at the last observed visit) and estimate their shared geographical pattern. The choice of the VL immunological outcome cut-off point was based on the global goal to achieve an undetectable VL to minimize viral transmission in the general population [15] while the choice of the CD4 cell counts immunological outcome cut-off point was based on earlier studies that have shown that individuals with CD4 cell counts of 200-350cells/uL immune deterioration patterns are not significantly different from those with CD4 cell counts less than 200cells/uL [5].

The MCAR model is a special type of Gaussian Markov random field prior models (GMRF) [14, 16]. Let the number of observed individuals with detectable VL (VL≥50copies/uL) beYVL, r, and the number of observed individuals with immune deterioration (CD4 < 350cells/uL) be YCD4, r at the last observed visit as reported for the region Sr, where the set of regions {Sr}, r = 1, 2, 3, …, R represents a finite number of the areas partitioned from the entire study area. In this study, we considered partitioning the Zimbabwe country into 10 provinces, i.e. R = 10. The geospatial variations are observed through an aerial map partitioned based on administrative spatial units. For each region, the expected number of individuals for each prognostic marker outcome (EVL, r and ECD4, r)was calculated using the following generic formula:

$$ {E}_{r,q}={n}_{r,q}{\overline{r}}_q\equiv {n}_{r,q}\left(\frac{\sum_r{y}_{r,q}}{\sum_r{n}_{r,q}}\right)\ \mathrm{for}\ \mathrm{r}=1,2,3,\dots, R\ \mathrm{and}\ \mathrm{q}= CD4, VL $$
(7)

where \( {\overline{r}}_q \) is the overall detectable VL (VL ≥ 50copies/uL) or immune deterioration (CD4 < 350 cell/uL) rate in the whole study region, nr, q is the at-risk population in the region r and yr, q is the total counts of individuals observed in the region. This approach is commonly referred to as “internal standardization”. The observed frequencies can be considered as the random variables while the expected frequencies are thought of as fixed and known functions of the at-risk population nr, q in the region r. In this study, we assumed that the observed count data follow a Poisson distribution for the two prognostic marker outcomes, i.e.

$$ {\displaystyle \begin{array}{l}{Y}_{q,r}\sim Poisson\left({E}_{q,r}{\theta}_{q,r}\right)\mathrm{for}\ \mathrm{r}=1,2,3,\dots, R\ \mathrm{and}\ \mathrm{q}= CD4, VL\\ {}\log \left({\theta}_{q,r}\right)=\log \left({E}_{q,r}\right)+{\alpha}_q+{S}_{q,r}\end{array}} $$
(8)

where αq is an intercept term representing the baseline (log) relative risk of disease progression outcomes q across the study region. We further assumed that the log relative risks are spatially correlated across the regions, and the log relative risks for the two prognostic markers are also correlated within each region r due to shared region-level unmeasured risk factors. These assumptions were supported through the intrinsic bivariate CAR prior of a 2 × R dimensional matrix of Sq, r values. The spatial prior is expressed as:

$$ {S}_q\mid {S}_{1\left(-q\right)},{S}_{2\left(-q\right)}\sim \mathrm{Bivariate}\ \mathrm{Normal}\left({\overline{\boldsymbol{S}}}_q,\boldsymbol{V}/{n}_q\right) $$
(9)

where S1(−q), S2(−q) denotes the elements of the 2 × R matrix, \( {\overline{\boldsymbol{S}}}_r=\left({\overline{S}}_{r,1},{\overline{S}}_{r,2}\right) \) and \( {\overline{S}}_{r,p}={\sum}_{\mathrm{f}\ \mathrm{in}\ {\sigma}_r}{S}_{fp}/{n}_r \) where σr and nr denote the set of labels of the “neighbours” of the region r and the number of neighbours, respectively assuming a p = 2. The matrix V is a 2 × 2 covariance matrix with diagonal elements v11 and v22 denoting the conditional variances of S1 and S2 respectively, and an off-diagonal element v12 denoting the conditional within-area covariance between S1 and S2.

The MCAR model was performed in the OpenBUGS, which is open-source statistical software and the code used is provided in the Additional file 1 section. We ran 10,000 Markov chain Monte Carlo (MCMC) simulations, burn-in of 2000 and thinning of 10. A prior sensitivity analysis was done through varying prior distributions and parameter values. We ran simultaneously two chains and the results reported are based on the better chain of initial values. Model diagnostics were assessed through trace plots which should traverse rapidly in the same region, density plots which should show a smooth curve and autocorrelation plots which should show a quick sharp drop in early lags.

Results

A total of 3896 participants contributed 8655 follow-up observations. There were 2551(65.5%) females and 1345(34.52%) were males. The average age was 38.23 ± 11.37 years and 1388(35.6%) participants were aged 35–44 years. At baseline, the median CD4 cell count was 211cells/uL with an interquartile range (IQR) of 114-320cells/uL and 1846 (47.4%) participants had CD4 cell counts of below 350cells/uL. The median VL was 57copies/uL (IQR: 48-66copies/uL) and 2557 (65.6%) participants had VL between 50 and 1000 copies/uL at baseline. All participants were on three-drug combination therapy and none were receiving protease inhibitors (PI) or integrase inhibitors INSTI) regimen. Most of the participants (95%) were on a first-line three-drug combination therapy which was a combination therapy of two nucleoside reverse transcriptase inhibitors (NRTI) namely Tenofovir (TDF) and Lamivudine (3TC); and one non-NRTI namely Efavirenz (EFV), that is, TDF + 3TC + EFV, while 5% were on second-line drugs.

Table 1 shows the simple regression model results used to generate the orthogonal covariates for both the CD4 cell counts and the VL model. We observed a significant correlation between CD4 cell counts and VL for both models and both slopes were negative as expected.

Table 1 Estimated regression coefficients for the simple linear regression model for the viral load on CD4 cell counts and CD4 cell counts on viral load

Table 2 presents a detailed summary of the log-likelihood values, LRT statistics, LRT p-values and the AIC values. Using the LRT results for the nested models, we assessed a better fit model between the no covariate model \( \left[{\lambda}_{ij}(t)={\lambda}_{ij(0)}\exp \left({\beta}_{ij}^{\prime}\right),i\ne j\right] \) and the orthogonal adjusted model [λij(t) = λij(0) exp(βij × P)] for both the CD4 cell counts and VL multistate models.

Table 2 Model selection process for the various model fitted for both viral load and CD4 cell counts multistate models

Both the CD4 cell count and the VL multistate models showed that the model with the orthogonal variable (model 2) was a better fit compared to the no variate model (model 1). Similarly, the multistate models which adjusted for the time-varying effects of the other prognostic marker, [λij, VL(t) = λij(0) exp(βij × CD4state)and λij, CD4(t) = λij(0) exp(βij × VLstate)] (model 3) were better fit compared to the no covariate models (model 1). Simultaneously adjusting for the orthogonal variable and the time-varying effects on both the CD4 cell count model and VL model (model 4) further improve the model, P-value< 0.05; hence, the interpretation of results was based on these models (model 4).

Modelling of CD4 cell counts adjusting for the viral load orthogonal and time-varying effects

The orthogonal and the time-varying VL effects were regressed in a single model. In this model, the movement between state i to state j for i > j defines an increase in CD4 cell counts which indicates immune recovery process while i < j defines a decrease in CD4 cell counts, which means immune deterioration process.

In Table 3, results show that the rates of decrease in CD4 cell counts (immune deterioration) among patients in state 2 (CD4 cell count range of 350–500 cells/uL; λ23 = 0.0404) were 11.22 times higher than the rates of increase in CD4 cell counts (immune recovery) among patients in state 3 (CD4 < 350cells/uL; λ32 = 0.0036). Patients in state 2 (CD4 cell count range of 350–500 cells/uL) were more likely to become LTFU (λ25 = 0.0551). The positive VL log-linear effects showed an increased risk of immune deterioration for an individual in state 2 (CD4 cell count range of 350–500 cells/uL) to state 3 (CD4 < 350cells/uL) \( \left({\beta}_{23}{P}_{VL(k)}^{\ast }=0.4762\right) \). Similarly, the log-linear effects of the time-varying VL variable showed an increased risk of immune deterioration for the λ23 transition (β23VLstates, k = 1.0791).

Table 3 The joint effects of time-varying viral load and orthogonal viral load variables on the CD4 cell counts transition rates Markov model

The time-varying VL values and the negative VL (residuals covariate) had an increasing effect on the risk of death in this cohort for individuals with CD4 < 350cell/uL (hazard ratio (HR) =1.75 and HR = 1.67, respectively). Regarding the time-varying effects of VL on the CD4 cell counts transition rates, as the CD4 cell counts increased for an individual with a CD4 cell count range of 350–500 cells/uL (λ23), the transition rates increased with an increase in VL levels over time. Similarly, individuals with CD4 ≥ 500cells/uL (λ15) and those with CD4 < 350cell/uL (λ35) showed an increased risk of becoming LTFU as their VL increased over time. A comparable trend was observed for individuals with CD4 < 350cell/uL (λ34) who exhibited an increased risk of death with an increase in VL over time.

Modelling of viral load adjusting for the CD4 cell counts orthogonal and time-varying effects

We fitted a time-homogeneous multistate Markov model to evaluate HIV disease progression based on VL defined states. The effects of the orthogonal and the time-varying CD4 cell counts effects were accounted for in a single model. In this model, the movement between state i to state j for i > j defines a VL suppression while i < j defines VL rebound.

In Table 4, the rates of viral rebound among patients with VL < 50copies/uL) (λ12 = 0.0286) were 1.78 times higher than the rates of VL suppression among patients with a VL range of 50-1000copies/uL (λ21 = 0.0161). However, patients with VL ≥ 1000copies/uL were 7.77 times more likely to die (λ32 = 0.0035 vs λ34 = 0.0272) and 9.8 times more likely to become LTFU (λ32 = 0.0035 vs λ35 = 0.0343). The time-varying CD4 cell counts had an increasing effect on the risk of viral increase from a VL range of 50-1000copies/uL to VL ≥ 1000copies/uL in this cohort, HR =1.23.

Table 4 The joint effects of time-varying CD4 cell counts levels and orthogonal CD4 cell counts variables on the viral load transition rates Markov model

Similarly, the negative CD4 cell count residual was associated with an increased risk of viral increase from a VL range of 50-1000copies/uL to VL ≥ 1000copies/uL \( \left({\beta}_{23}{P}_{CD4(k)}^{\ast }=0.6052;\mathrm{HR}=1.83\right) \); and an increased risk of death among individuals with VL ≥ 1000copies/uL \( \left({\beta}_{34}{P}_{CD4(k)}^{\ast }=0.3475;\mathrm{HR}=1.42\right) \). Regarding the time-varying effects of CD4 cell counts on the VL transition rates, as the VL increased (λ12 and λ23) the transition rates increased with a decrease in the CD4 cell counts over time while as the VL decreased (λ21, λ31 and λ32) the transition rates decreased with a decrease in the CD4 cell counts over time. The mortality rates of individuals with VL < 50copies/uL and those with VL range between 50-1000copies/uL increased over time as CD4 cell count decreases (λ14 and λ24). Regardless of the individual’s VL state, the transition rates to become LTFU decreased with a decrease in CD4 cell counts. Individuals with high CD4 cell counts were more likely to become LTFU in this cohort (λ15, λ25 and λ35); however, patients with CD4 ≥ 500cells/uL and a VL range of 50-1000copies/uL had an increased risk of becoming LTFU.

Multistate Markov models assessment

We performed a post-estimation test to assess which of the two multistate models fit better on this data in describing the HIV disease progression, mortality and becoming LTFU. We used the prevalence plot to compare the observed and the expected percentage prevalence for CD4 cell counts and VL multistate models. The VL multistate model showed a perfect fit for state 3(VL ≥ 1000copies/uL) and fair (moderate) fit for state 1 (VL < 50copies/uL) and state 4 (died). State 3 (VL ≥ 1000copies/uL) and state 5 (LTFU) of this model were not perfectly fitted by this model (Fig. 2).

Fig. 2
figure 2

The percentage prevalence plot of orthogonal CD4 cell counts and the time-varying CD4 cell counts values on the viral load multistate model  [State 1(VL < 50), State 2 (50 ≤ VL < 1000), State 3(VL ≥ 1000), State 4 (dead) and State 5 (Loss to follow-up (LTFU)]

Regarding the CD4 cell counts multistate model, state 1 (CD4 ≥ 500cell/uL) was perfectly fitted while state 2 (350 ≤ CD4 < 500) and state 4 (died) were fairly fitted. State 3 (CD4 < 350cell/uL) was perfectly fitted between time 0 and time 5 years only while state 5 (LTFU) was not perfectly fitted (Fig. 3). Since the percentage prevalence plots could not give conclusive results on the better model and the LRT was not appropriate as the models were not nested, we used the AIC values. The VL multistate model had an AIC = 7002.71 while the CD4 cell counts multistate model had an AIC = 7337.46; hence, the VL multistate model was a better fit model.

Fig. 3
figure 3

The percentage prevalence plot of orthogonal viral load variable and the time-varying viral load variables on the CD4 multistate model [State 1(CD4 ≥ 500), State 2 (350 ≤ CD4 < 500), State 3 (CD4 < 350), State 4 (dead) and State 5 (Loss to follow-up (LTFU)]

Sub-analysis of the spatial covariate (province of ART enrolment) on the multistate models and joint mapping of the two immunological HIV conditions using the multivariate intrinsic CAR model

We further assessed the effects of the region variable (province) on the VL and CD4 cell count multistate models to identify any association between the observed transition processes and the region. We estimated the HR summarised in Table 5. In general, we observed that there was an association between the transition intensities and the province administrative level. We found that 57.6% of the transition rates for CD4 cell count and VL models adjusting for province covariate showed the same effect direction (risk or protective effect) of the hazard ratios; however, there were variations in the magnitude of the risk of transitions among provinces. Immune deterioration (decrease in CD4 cell counts) was evident among patients belonging to Mashonaland East (P4) (state 1 to 3 and state 2 to 3) and Matabeleland North (P7) (state 1 to 2), while the immune recovery was most evident among patients from Masvingo (P6).

Table 5 The spatial effects (province-level) on the viral load and CD4 cell counts multistate model adjusted for time-varying effects and residual effect of the other prognostic marker with Harare province as a reference category

Becoming LTFU was observed to be highest among patients from Mashonaland East (P4) regardless of their CD4 cell count. We found that the risk of viral rebound to above 50copiels/uL was high among patients from Matabeleland North (P7) (state 1 to 2), Mashonaland East (P4) (state 1 to 3) regions. In contrast, viral suppression to undetectable levels (VL < 50copies/uL) was evident among patients from Masvingo (P6) (state 2 to 1) and Matabeleland North (P7) (state 3 to 1). The risk of becoming LTFU was high among patients from Mashonaland East (P4) if their VL < 50copiels/uL, Masvingo (P6) if their VL range was between 50-1000copies/uL and Mashonaland West (P5) if their VL ≥ 1000copies/uL.

To get a pictorial view of the spatial patterns and correlation between the CD4 cell counts marker and the VL measurements, we fitted the multivariate intrinsic CAR prior model with the province as the spatial unit. We jointly modelled those patients who had a VL ≥ 50copies /uL (VL state 2 and 3 combined) at the end of the follow-up to define that group that might not have attained undetectable VL or have a VL rebound to detectable levels, and those patients who had a CD4 < 350cells/uL (CD4 state 3) to define that group that is still in the immune-deterioration phase at the end of the follow-up period.

Table 6 shows the posterior estimates after the joint mapping of the two immunological outcomes for HIV disease progression among ART patients based on CD4 cell counts and VL. The posterior correlation between the spatially structured risk components of having a detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) was 61.3% (95% credible interval (CI): 47–97%). This strong correlation suggests strong shared geographical patterns of the risk of immune deterioration defined by these two prognostic markers. The baseline (log) relative risk of the HIV disease progression based on detectable VL (VL ≥ 50copies/uL) was estimated at − 0.472 while that of immune deterioration (CD4 < 350cells/uL) was − 0.043.

Table 6 The posterior estimates for the multivariate conditional autoregressive model for the two HIV disease progression prognostic markers

The joint mapping of the posterior relative risk of the HIV disease progression defined by the VL (RR1) and CD4 cell counts (RR2) under the MCAR model is shown in Fig. 4. The dark blue or deep grey colours indicate areas with high relative risk while the light blue or light grey colours show regions with the lowest relative risks of the defined conditions. These maps show a geographical overlap of detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) in seven provinces. Patients in Matabeleland North province bordering Botswana and Zambia; and Mashonaland East province toward the Mozambique border had higher relative risks of detectable VL (VL ≥ 50copies/uL), and immune deterioration (CD4 < 350cells/uL) outcomes in addition to Bulawayo and Harare metropolitan provinces.

Fig. 4
figure 4

The joint mapping of the posterior relative risk of the HIV disease progression defined by the viral load (RR1) and CD4 cell counts (RR2) under the multivariate conditional autoregressive model. The maps were generated from OPENBUGS version 3.2.3 [17] https://www.mrc-bsu.cam.ac.uk/software/bugs/openbugs/

Discussion

The main objective of our study was to jointly model HIV disease progression using two 5-stage time-homogenous multistate Markov models based on CD4 cell counts and VL. The multistate models accounted for the orthogonal and the time-varying covariates simultaneously which has been a gap in previous studies. The spatial overlap of the HIV disease progression outcomes (detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) at the last observed visit conditional not to have died or become loss to follow-up (LTFU) was also described.

We fitted different time-homogeneous multistate Markov models and observed that multistate models with both the orthogonal and the time-varying variables were the best fitting models to describe the HIV disease progression. The VL multistate model was superior in predicting the HIV disease progression patterns in this cohort compared to the CD4 cell count multistate model. These findings support what earlier studies have reported on the superiority of VL in monitoring HIV disease progression patterns compared to CD4 cell counts [7, 18]. Our findings further support the global guidelines on the use of VL as the primary laboratory marker to routinely monitor HIV disease progression amongst ART patients.

The VL multistate model with CD4 cell counts in the covariate matrix showed that VL rebound and VL increase transition rates (state 1 to 2; state 1 to 3 and state 2 to 3) were higher than the VL suppression transition rates (state 2 to 1; state 3 to 2 and state 3 to 1). This means viral rebound was more likely compared to viral suppression in this cohort. We also found that as the VL increased (λ12 and λ23) the transition rates increased with a decrease in the CD4 cell counts over time while as the VL decreased (λ21, λ31 and λ32) the transition rates decreased with a decrease in the CD4 cell counts decreased over time. This finding aligns with the negative correlation which exists between these two prognostic markers that as CD4 cell count levels decreases, VL increases. These finding could be explained by poor ART adherence which remains a challenge among ART patients [19]. Moreover, treatment failure or treatment side-effects may be possible underlying factors to explain these results [18]. With effective and potent ART; and without non-adherence challenges, the VL is expected to decrease as the CD4 cell counts increase. Therefore, interventions that support ART adherence like support groups, quick identification of treatment failure and continual education of the effectiveness of ART should be strengthened.

Becoming LTFU and mortality were more likely among patients with a VL ≥ 1000copies/uL (state 3 to 4 and state 3 to 5). Patients with CD4 cell counts lower than expected (negative orthogonal variable) were associated with an increased VL and mortality if VL > 50copies/uL. This means as patients with high VL immune deteriorate, the risk of death becomes higher. These findings support other earlier studies that once the immune system deteriorates, the chances of immune recovery become slim, and thus when most deaths occur [20,21,22]. The risk of becoming LTFU was high among patients with VL ≥ 50copies/uL and CD4 < 350cells/uL. This can be explained by the fact that more ill patients tend to drop out from the ART programme and are classified as LTFU [5]. However, some of these patients if tracked, would have died, which may result in a misclassification error of the LTFU outcome. In contrast, patients with higher CD4 cell counts were associated with becoming LTFU. This could be explained by the fact that much sicker patients are more likely to be bedridden; hence, they are hospital-bound while healthier patients may become LTFU as a result of “silent-transfers” to nearby health [23, 24].

Adjusting for the spatial covariate in the VL multistate model, as expected, the regions with a high risk of viral rebound were also more likely to have a high risk of low CD4 cell counts. The intrinsic MCAR prior model further confirmed a strong overlapping geographical correlation between individuals with a detectable VL (VL ≥ 50copies/uL) and immune deterioration (CD4 < 350cells/uL) of 66%; hence, a shared geographical pattern of relative risks between these two outcomes exists. Patients staying in provinces that border with nearby countries had high relative risks of immune deterioration and detectable VL, particularly, those from Matabeleland North province in the northern part of Zimbabwe. Matabeleland North province has a busy truck route to neighbouring countries and high mobility of local and tourists. Earlier studies looking at spatial heterogeneity of viral suppression in this province reported similar results [8, 25]. People in this province are likely to present late for health care [26], delay ART initiation [27] or engage in sexual activities with multiple partners which subsequently compromise viral suppression among HIV patients [28, 29]. Therefore, interventions such as self-testing [30], pre-exposure prophylaxis (PrEP) pills for high-risk groups [31] should be intensified in such regions.

The reported results should be inferred in light of some limitations. Firstly, the dataset used in this study is very small and might not be an accurate representation of the Zimbabwean HIV population. We included only health facility data liked to the ePMS with both VL and CD4 cell counts repeated measurements. In this regard, the MCAR assumption might not have been fully satisfied. The spatial effects were observed at a higher level which might not precisely pinpoint the marginalized areas to guide policy. This was as a result of our small sample size in this study; however, future studies should consider lower administrative level spatial units. Moreover, the MCAR model only describes the spatial interaction across the error terms to explain spatial autocorrelation; however, this model falls short when a more direct presentation of spatial interaction is desired. The VL and CD4 cell counts measurements were not randomly done rather differential monitoring was implemented due to resource constraints. The joint model could not account for the biological order of the association between VL and CD4 cell count, that is, VL change may precede CD4 change. We could not adjust for ART adherence, comorbidities (tuberculosis, diabetes and hypertension); and demographic characteristics (age and gender) due to model convergence issues. Multistate Markov models can estimate multiple transition rates and outcomes simultaneously compared to the Cox proportional hazard models. However, the assumption of constant hazard function does not reflect reality and also the Markov process has a memory loss property which may be a limitation in HIV studies [6]. Despite these limitations, this study managed to provide useful information in HIV monitoring through jointly modelling two HIV prognostic markers and identifying regions with poor immunological outcomes.

Conclusion

In conclusion, the findings from this study provide a foundation on the HIV disease progression patterns in Zimbabwe to guide policy going forward and motivate future statistical modelling to consider such modelling approaches in infectious disease to guide policy and programme management. We found that VL models were much more superior compared to CD4 cell count models in HIV monitoring. We also observed that after joint modelling of CD4 cell counts and VL in a single model, the rates of VL rebound risks are still higher than viral suppression. We observed that having a high VL (VL ≥ 1000copies/uL) increased the risk of becoming LTFU and death. The spatial overlap between VL and CD4 cell counts indicates the inter-linkages between these two markers in HIV monitoring and the shared geographical correlation of HIV disease progression immunological outcomes was high. With the global efforts to achieve zero HIV incidence, interventions which encourage ART adherence like support groups and continual educational support on the barriers for ART uptake is crucial. Region-specific interventions which focus on behavioural change, social support, increased uptake of HIV preventions like PrEP should be strengthened among high-risk groups all in the quest to achieve zero HIV incidences by 2030.

Availability of data and materials

The data that supports the findings of this study are available through an application process from the Ministry of Health and Child Care in Zimbabwe which is the custodian of the ePMS data through the AIDS/TB Unit who manages and oversees the ePMS data collection process and are the custodian of the database.

References

  1. World Health Organization. Consolidated guidelines on the use of antiretroviral drugs for treating and preventing HIV infection: recommendations for a public health approach: World Health Organization; 2016. p. 1–155. Available from: https://www.who.int/hiv/pub/arv/arv-2016/en/ (Accessed 15 Jan 2020)

    Google Scholar 

  2. World Health Organization. HIV Treatment and Care What’s new in Treatment Monitoring: Viral load and CD4 Testing: WHO; 2017. p. 1–2. Available from: https://www.who.int/hiv/pub/arv/treatment-monitoring-info-2017/en/ Accessed 26 November 2018

    Google Scholar 

  3. Meira-Machado L, De Uña-álvarez J, Cadarso-suárez C, Andersen PK. Multi-state models for the analysis of time-to-event data. Stat Methods Med Res. 2009;18(2):195–222. https://doi.org/10.1177/0962280208092301.

    Article  PubMed  Google Scholar 

  4. Matsena Zingoni Z, Chirwa TF, Todd J, Musenge E. A review of multistate modelling approaches in monitoring disease progression: Bayesian estimation using the Kolmogorov-chapman forward equations. Stat Methods Med Res. 2021;962280221997507(5):1373–92. https://doi.org/10.1177/0962280221997507.

    Article  Google Scholar 

  5. Matsena Zingoni Z, Chirwa TF, Todd J, Musenge E. HIV disease progression among antiretroviral therapy patients in Zimbabwe : a multistate Markov model. Front Public Health. 2019;7:1–15. https://doi.org/10.3389/fpubh.2019.00326.

    Article  Google Scholar 

  6. Shoko C, Chikobvu D. Time-homogeneous Markov process for HIV / AIDS progression under a combination treatment therapy: cohort study, South Africa. BMC Theor Biol Med Model. 2018;15(1):1–14. https://doi.org/10.1186/s12976-017-0075-4.

    Article  Google Scholar 

  7. Shoko C, Chikobvu D. A superiority of viral load over CD4 cell count when predicting mortality in HIV patients on therapy. BMC Infect Dis. 2019;19(1):1–10. https://doi.org/10.1186/s12879-019-3781-1.

    Article  Google Scholar 

  8. Matsena Zingoni Z, Chirwa TF, Todd J, Musenge E. Spatial Heterogeneity of Viral Suppression and Viral Rebound Patterns among ART Patients in Zimbabwe from 2004 to 2017 : A Bayesian Mixed Effects Multistate Model. Int J Stat Med Res. 2019;8:98–113. https://doi.org/10.6000/1929-6029.2019.08.13.

    Article  Google Scholar 

  9. Ford N, Roberts T, Calmy A. Viral load monitoring in resource-limited settings. AIDS. 2012;26(13):1719–20. https://doi.org/10.1097/QAD.0b013e3283543e2c.

    Article  PubMed  Google Scholar 

  10. Odhiambo C, Kareko MJ. An evaluation of frequentist and Bayesian approach to geo-spatial analysis of HIV viral load suppression data. Int J Stat Appl. 2019;9(6):171–9. https://doi.org/10.5923/j.statistics.20190906.01.

    Article  Google Scholar 

  11. Shoko C, Chikobvu D, Bessong PO. A Markov model to estimate mortality due to HIV/AIDS using viral load levels-based states and CD4 cell counts: a principal component analysis approach. Infect Dis Ther. 2018;7(4):457–71. https://doi.org/10.1007/s40121-018-0217-y.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Ministry of Health and Child Care. Electronic Patient Management System (ePMS)- Zimbabwe. 2012.

    Google Scholar 

  13. Mukumba T. Usability assessment of the electronic patient management system for AIDS & TB Services in Zimbabwe. Open Access Library J. 2014;1(08):1–18. https://doi.org/10.4236/oalib.1101119.

    Article  Google Scholar 

  14. Besag J, York J, Mollie A. Bayesian image restoration, with two applications in spatial statistics. Ann Inst Stat Math. 1991;43(1):1–59. https://doi.org/10.1007/BF00116466.

    Article  Google Scholar 

  15. UNAIDS. Undetectable = Untransmittable Public Health and HIV viral load suppression; 2018. p. 1–4. Available from: http://www.unaids.org/sites/default/files/media_asset/undetectable-untransmittable_en.pdf Accessed 5 Feb 2019

    Google Scholar 

  16. Leroux BG, Lei X, Breslow N. Estimation of Disease Rates in Small Areas: A new Mixed Model for Spatial Dependence. Statistical Models in Epidemiology, the Environment, and Clinical Trials. 2000;(1999):179–191. doi: https://doi.org/10.1007/978-1-4612-1284-3_4

    Book  Google Scholar 

  17. Thomas A. OpenBUGS: Constructing Mcmc Software: John Wiley & Sons Incorporated; 2007. Wiley Series in Computational Statistics Series

    Google Scholar 

  18. Nicholas S, Poulet E, Wolters L, Wapling J, Rakesh A, Amoros I, et al. Point-of-care viral load monitoring: outcomes from a decentralized HIV programme in Malawi. J Int AIDS Soc. 2019;22(8):1–9. https://doi.org/10.1002/jia2.25387.

    Article  Google Scholar 

  19. Rachlis B, Burchell AN, Gardner S, Light L, Raboud J, Antoniou T, et al. Social determinants of health and retention in HIV care in a clinical cohort in Ontario, Canada. AIDS Care. 2017;29(7):828–37. https://doi.org/10.1080/09540121.2016.1271389.

    Article  PubMed  Google Scholar 

  20. Grover G, Gadpayle AK, Swain PK, Deka B. A Multistate Markov Model Based on CD4 Cell Count for HIV/AIDS Patients on Antiretroviral Therapy (ART). Int J Stat Med Res. 2013;2(1998):144–51. https://doi.org/10.6000/1929-6029.2013.02.02.08.

    Article  Google Scholar 

  21. Dessie ZG. Multistate models of HIV / AIDS by homogeneous semi-Markov process. Am J Biostat. 2014;4(2):21–8. https://doi.org/10.3844/ajbssp.2014.21.28.

    Article  Google Scholar 

  22. Goshu AT, Dessie ZG. Modelling progression of HIV/AIDS disease stages using semi-Markov processes. J Data Sci. 2013;11(2):269–80.

    Article  Google Scholar 

  23. Geng EH, Odeny TA, Lyamuya R, Nakiwogga-Muwanga A, Diero L, Bwana M, et al. Retention in care and patient-reported reasons for undocumented transfer or stopping care among HIV-infected patients on antiretroviral therapy in eastern Africa: application of a sampling-based approach. Clin Infect Dis. 2016;62(7):935–44. https://doi.org/10.1093/cid/civ1004.

    Article  PubMed  Google Scholar 

  24. Mutasa-Apollo T, Shiraishi RW, Takarinda KC, Dzangare J, Mugurungi O, Murungu J, et al. Patient retention, clinical outcomes and attrition-associated factors of HIV-infected patients enrolled in Zimbabwe’s National Antiretroviral Therapy Programme, 2007-2010. PLoS One. 2014;9(1):2007–10. https://doi.org/10.1371/journal.pone.0086305.

    Article  CAS  Google Scholar 

  25. Zimbabwe Ministry of Health and Child Care (MOHCC). Zimbabwe Population-Based HIV Impact Assessment (ZIMPHIA) 2015–2016: First Report. MOHCC Harare; 2017. p. 1–79. Available from: https://phia.icap.columbia.edu/wp-content/uploads/2017/11/ZIMPHIA_First_Report_FINAL.pdf Accessed 7 Dec 2019

    Google Scholar 

  26. Nyika H, Mugurungi O, Shambira G, Gombe NT, Bangure D, Mungati M, et al. Factors associated with late presentation for HIV / AIDS care in Harare City, Zimbabwe, 2015. BMC Public Health. 2016;16(369):1–7. https://doi.org/10.1186/s12889-016-3044-7.

    Article  CAS  Google Scholar 

  27. Maponga BA, Chirundu D, Gombe NT, Tshimanga M, Bangure D, Takundwa L. Delayed initiation of antiretroviral therapy in TB/HIV co-infected patients, Sanyati district, Zimbabwe, 2011–2012. Pan Afr Med J. 2015;21:2011–2. doi: https://doi.org/10.11604/pamj.2015.21.28.5195

  28. Wabiri N, Shisana O, Zuma K, Freeman J. Assessing the spatial nonstationarity in relationship between local patterns of HIV infections and the covariates in South Africa: a geographically weighted regression analysis. Spatial Spatio-temporal Epidemiol. 2016;16:88–99. https://doi.org/10.1016/j.sste.2015.12.003.

    Article  Google Scholar 

  29. Rangarajan S, Donn JC, Giang LT, Bui DD, Hung Nguyen H, Tou PB, et al. Factors associated with HIV viral load suppression on antiretroviral therapy in Vietnam. J Virus Erad. 2016;2(2):94–101. https://doi.org/10.1016/S2055-6640(20)30466-0.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Sharma M, Ying R, Tarr G, Barnabas R. Systematic review and meta-analysis of community and facility-based HIV testing to address linkage to care gaps in sub-Saharan Africa. Nature. 2015;528(7580):S77–85. https://doi.org/10.1038/nature16044.

    Article  PubMed  PubMed Central  Google Scholar 

  31. World Health Organization. Guideline on When To Start Antiretroviral Therapy and on Pre-Exposure Prophylaxis for HIV. Geneva: World Health Organization; 2015. p. 1–78. Available from: http://www.who.int/hiv/pub/guidelines/earlyrelease-arv/en/. Accessed 5 May 2018

Download references

Acknowledgements

Our acknowledgements go to the Ministry of Health and Child Care, AIDS/TB Units department for the support and compilation of the data used in this study. We also thank the Division of Epidemiology and Biostatistics at the School of Public Health for their support. Thanks to Dr. Claris Shoko for technical assistance on the data analysis stage.

Funding

This study was supported by the Developing Excellence in Leadership, Training and Science (DELTAS) Africa Initiative Sub-Saharan Africa Consortium for Advanced Biostatistics (SSACAB [Grant No. DEL-15-005]. The DELTAS Africa Initiative is an independent funding scheme of the African Academy of Sciences (AAS) Alliance for Accelerating Excellence in Science in Africa (AESA) and is supported by the New Partnership for Africa’s Development Planning and Coordinating Agency (NEPAD Agency) with funding from the Wellcome Trust [Grant No. 107754/Z/15/Z] and the United Kingdom government. The views expressed in this publication are those of the authors and not necessarily those of the AAS, NEPAD Agency, Wellcome Trust, the UK government or the Ministry of Health and Child Care, Zimbabwe since these parties did not play any role in the design of the study, data analysis, interpretation of results and manuscript writing. This work was also supported by the Wits University Research Committee (URC) [Grant number: URC-2021].

Author information

Authors and Affiliations

Authors

Contributions

ZMZ cleaned and analyzed the data and drafted the manuscript, JT and TFC reviewed the manuscript and advised on analysis, EM guided and oversaw the analysis and reviewed the manuscript. All authors read and approved the final manuscript for submission.

Corresponding author

Correspondence to Zvifadzo Matsena Zingoni.

Ethics declarations

Ethics approval and consent to participate

The permission to use the routinely collected HIV data amongst ART patients was sort and granted from the Zimbabwe Ministry of Health and Child Care through the AIDS/TB Unit who oversees the data collection and storage. Ethical approval was granted by the University of Witwatersrand’s Human Research Ethics Committee (Medical) (Clearance Certificate No. M170673).

Consent for publication

Not applicable. All materials in this manuscript are mine and were not extracted from any existing sources.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Matsena Zingoni, Z., Chirwa, T.F., Todd, J. et al. Markov modelling of viral load adjusting for CD4 orthogonal variable and multivariate conditional autoregressive mapping of the HIV immunological outcomes among ART patients in Zimbabwe. Theor Biol Med Model 18, 16 (2021). https://doi.org/10.1186/s12976-021-00145-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12976-021-00145-y

Keywords