Differences in Transcriptional Dynamics Between T-cells and Macrophages as Determined by a Three-State Mathematical Model

HIV-1 viral transcription persists in patients despite antiretroviral treatment, potentially due to intermittent HIV-1 LTR activation. While several mathematical models have been explored in the context of LTR-protein interactions, in this work for the first time HIV-1 LTR model featuring repressed, intermediate, and activated LTR states is integrated with generation of long (env) and short (TAR) RNAs and proteins (Tat, Pr55, and p24) in T-cells and macrophages using both cell lines and infected primary cells. This type of extended modeling framework allows us to compare and contrast behavior of these two cell types. We demonstrate that they exhibit unique LTR dynamics, which ultimately results in differences in the magnitude of viral products generated. One of the distinctive features of this work is that it relies on experimental data in reaction rate computations. Two RNA transcription rates from the activated promoter states are fit by comparison of experimental data to model predictions. Fitting to the data also provides estimates for the degradation/exit rates for long and short viral RNA. Our experimentally generated data is in reasonable agreement for the T-cell as well macrophage population and gives strong evidence in support of using the proposed integrated modeling paradigm. Sensitivity analysis performed using Latin hypercube sampling method confirms robustness of the model with respect to small parameter perturbations. Finally, incorporation of a transcription inhibitor (F07#13) into the governing equations demonstrates how the model can be used to assess drug efficacy. Collectively, our model indicates transcriptional differences between latently HIV-1 infected T-cells and macrophages and provides a novel platform to study various transcriptional dynamics leading to latency or activation in numerous cell types and physiological conditions.

of Histone H1 [41][42][43] . Here, we immunoprecipitated CDK9/Cyclin T 1 complexes from infected T-cell whole cell extracts and used samples from 0 and 6 hrs (+/− inducer) to better represent the transitions of the HIV-1 LTR with the absence of an inducer representing the transition to an intermediate state (LTR R to LTR I ) and the presence of an inducer representing the transition to an activated state (LTR I to LTR A ). We used Histone H 1 as a substrate as it has previously been shown to serve as a substrate for the p-TEFb complex 36,37 . Results in Fig. 2b show a background level of CDK9 kinase activity at 0 h (lane 1), which may indicate presence of short abortive transcripts, followed by a 34.7% increase at 6 h (lane 2). Similar results were obtained when using a whole cell extract from induced HIV-1 infected cells (20.5%; lane 3). Collectively, these data indicate that the transient transcriptional activation in HIV-1 infected T-cells that produces short HIV-1 RNA transcripts and occasional genomic readthroughs may potentially be a more rapid event allowing the LTR to return to a more quiescent state quickly as compared to the sustained Tat-activated transcription that results in production of infectious virus.
To construct parameter values for the mathematical model ( Fig. 1), densitometry counts of the kinase assay ( Fig. 2a,b) were taken to determine the rate of transition between each of the HIV-1 LTR states (LTR R , LTR I , LTR A ). Briefly, to determine the parameter values k ON (LTR R to LTR I ) and k A Tat (LTR I to LTR A ), we utilized densitometry counts of the CDK9-IP kinase assays between 0 and 6 h (Fig. 2b). The transition from LTR R to LTR I is measured in the absence of an inducer, while the transition from LTR I to LTR A is measured in the presence of an inducer of viral transcription. Furthermore, the reverse rates (k OFF and k I ; Fig. 1) utilize densitometry counts of the total kinase assay (Fig. 2a) between 6 and 24 h in the absence and presence of an inducer, respectively. Densitometry analysis was normalized as described in Materials and Methods Section.
The results in Fig. 2c demonstrate the proportion of T-cell HIV-1 LTRs which are in a repressed state at a given time over 120 h, as predicted by our mathematical model. The black line demonstrates the value of the original parameter set, while the grey lines are all the realizations with respect to the sampling of parameters using a Latin hypercube sampling method. The dashed green, red and blue lines represent 80%, 90% and 95% confidence intervals, respectively. These results show a sharp decline in the proportion of HIV-1 LTRs in the repressed state following activation, eventually stabilizing at low proportions (1.99%), as expected. Similarly, the results in  The HIV-1 LTR is represented in three states of activation (as opposed to a two-stage model where activation by noise and its regulation in Tat activity is not considered); LTR R denotes a repressed state (i.e. latency); LTR I represents an intermediate state of activation; and LTR A is a Tat-dependent activated state of the HIV-1 LTR in which full viral production is possible. The terms k ON and k OFF represent the rate of activation from latency and the return to latency, respectively. K A (Tat) is the rate of Tat-dependent activation from an intermediate LTR to a fully activated LTR state and the term k I represents the rate in the opposite direction. The diagram depicts the creation of two species of HIV-1 RNAs termed TAR and env (envelope). The rate at which TAR RNA is created is given by α m1R , α m1I and, α m1A and the TAR degradation/exportation rate is denoted by γ m1 . The HIV-1 RNA species env (genomic) is produced by the intermediate state LTR (env I ) at a rate of α m2I . HIV-1 Tat results in a positive feedback mechanism for virus replication, production at a rate of α p1 and Pr55 (Gag) production at a rate of α p2 . The activated state of the LTR also produces envelope, termed env A , at a Tat-dependent rate of f m2 (Tat). The Tat-dependent production of env A also results in the production of Pr55 (Gag) at a rate α p2 . Finally, this model depicts the rate cleavage of Pr55 (Gag) to the viral protein p24 via α p3 . Degradation rates (γ p1 , γ p2 , γ m1 , γ m2 ) were set to zero for preliminary results.
likely responsible for the persistent transcription of HIV-1 RNAs seen in long-term, cART treated patients 6,32,44 . Interestingly, despite vastly different approaches, these findings are in line with a model described by Razooky et al. which encompasses a two-state model of the HIV-1 LTR that incorporates Tat as a feedback mechanism. Both the model presented here and that which is proposed by Razooky and colleagues show a transient state, or a state with low probability, of LTR activation that is independent of Tat 24,45 . The proportion of HIV-1 LTRs in the activated state (LTR A ) is shown in Fig. 2e. As expected, the proportion of LTR A steadily increased after activation with 20% FBS media and treatment with an inducer, which resulted in the production of full length, genomic HIV-1 RNA and the production of infectious virions with approximately 92.64% of LTRs in an active state at 120 h. Collectively, the relative proportions of LTR activation states changes over time after induction of the virus in T-cells is presented in Fig. 2f, indicating that there is dynamic activity over a range of time with the most diversity of LTR states occurring between approximately 20-24 h.

RNA production in T-cells as a function of LTR activation.
We next examined the production of two HIV-1 viral transcripts including TAR RNA, a short, non-coding viral RNA produced as a result of non-processive transcription 46 , and mRNA responsible for the production of the HIV-1 envelope protein, termed env, as a marker of full length genomic transcription. As previously described, 3 biological replicates of infected T-cells were cultured in low serum media (0.1% FBS) for 36 h, and then were placed in 20% FBS media. Samples were analyzed at 0 and 24 h following the addition of high serum media by RT-qPCR in technical triplicate for the production of TAR and env RNA to determine the production parameters from LTR I . Figure 3a shows a high level of both TAR and env RNA production associated with a repressed transcriptional state (LTR R , 0 h), with an  1 infected T-cells) were cultured in low serum media (0.1% FBS) for 36 h, and subsequently incubated in 20% FBS media and treated with an inducer (IR). Resulting whole cell extracts were analyzed using an in vitro total kinase assay (a) or a CDK9 IP kinase assay (b) to assess for changes in the HIV-1 LTR. Biochemical data was used to construct parameters for mathematical modeling to determine relative proportions of the HIV-1 LTR in the various states; repressed (c), intermediate (d), and activated (e) over 120 h. The black line demonstrates the solved value of the original parameter set, while the grey lines are all the realizations with respect to the sampling of parameters using a Latin hypercube sampling method. The dashed green, red and blue lines represent 80%, 90% and 95% confidence intervals, respectively. (f) Overlay of all three LTR states; repressed (LTR R , red), intermediate (LTR I , blue), activated (LTR A , green). Selected lanes from the same blot with identical exposure settings are presented in panel a. were cultured in 0.1% serum media for 36 h, followed by incubation in 20% FBS media and treated with an inducer (PMA/PHA) at 0 h in biological triplicate. Samples were analyzed by RT-qPCR in technical triplicate at 0, 24, 48, 72, 96, and 120 h post-induction for TAR and genomic RNA. Values were used to fit parameter rates for TAR (c) and env (d) production from LTR A . Obtained RT-qPCR values are shown in conjunction with model predictions using solved fitted parameters. Parameter rates for mathematical modeling were constructed using RT-qPCR values and numerical estimations. Constructed parameters were used to mathematically model production of TAR RNA (e), env RNA from LTR I (f), Tat protein (g), and env RNA from LTR A (h). The black line represents the solved value of the original parameter set, the grey lines are all the realizations with respect to the sampling of parameters using a Latin hypercube sampling method. The dashed green, red and blue lines represent 80%, 90% and 95% confidence intervals, respectively. www.nature.com/scientificreports www.nature.com/scientificreports/ average of 3.00 × 10 6 and 8.31 × 10 4 , respectively, from three biological replicates analyzed in technical triplicate. To accurately estimate RNA production from LTR A , given that the modeled time frame is 120 h, T-cells were cultured as described above with induction using PMA/PHA at time 0. Samples analyzed by RT-qPCR for viral RNA production at 0, 24, 48, 72, 96, and 120 h to determine parameter values for the production of RNA from LTR I and LTR A . Similar to Fig. 3a, the results in Fig. 3b also show a high TAR and env RNA associated with a repressed transcriptional state (LTR R , 0 h), consisting of an average of 1.45 × 10 7 and 7.08 × 10 6 , respectively, from three biological replicates analyzed in technical triplicate. Taken together, these results suggest that there is an increase in transcription of both short, noncoding RNAs and long genomic RNAs associated with the LTR activation.
To model experimental results mathematically, we determined RNA production rates for each of the defined RNA parameters (Fig. 1). In our model the repressed state of the LTR (LTR R ) does not result in the production of infectious virus as demonstrated by no transcription of full-length genomic RNA (env). However, there is persistent, non-processive transcription which results in the production of TAR RNA, as supported by numerous studies which show the presence of TAR RNA in both in vitro and in vivo samples despite suppressive antiretroviral regimens 6,44,47,48 . These studies have found that intracellular levels during antiretroviral treatment can range from 1 × 10 1 to 1 × 10 5 TAR RNA copies in vivo 6,44 . Along these lines, we have estimated the production of TAR RNA in infected T-cells to occur at a rate of 2.5 × 10 4 copies/mL/h (α m1R ). The hourly rates at LTR I for TAR RNA production (α m1I ) and env RNA production (α m2I ) were constructed using experimentally derived averages at the corresponding time ( Table 1). The parameter values for TAR and env RNA production from LTR A , α m1A and α m2A , (Fig. 3c,d) as well as degradation/exit rates of both RNAs were fit to the experimental data in Fig. 3b using the standard least square fitting procedure (see Materials and Methods). The degradation/exit rates can encompass multiple potential mechanisms beyond just degradation. This can include the export of viral RNA out of the cell through such cellular secretory pathways as secretory autophagy and extracellular vesicle release 33,35,44,48,49 . Furthermore, to correlate RNA transcription rates with the expression of HIV-1 Tat protein, we utilized densitometry counts of Tat Western blots (Supp. Fig. 1) at the corresponding times to calculate a Tat production parameter value (α p1 ) which is modeled as a function of env RNA production from LTRs in the intermediate state. The production of HIV-1 Tat from the repressed and activated LTR state is not included in our model due to predominant production of Tat from early doubly spliced mRNAs which are produced from full length transcription of genomic RNA. Detailed information regarding governing equations is presented in the Materials and Methods section.
According to our mathematical model, the production of TAR RNA is low during the repressed state but upon progression through the intermediate to the active state of the LTR there was sustained, increasing production of TAR RNA (Fig. 3e). The dynamics of env (genomic) RNA production from LTR I (Fig. 3f) show a slight increase in production starting at approximately 1.49 h and ending at approximately 81.18 h, which correlates to the modeled LTR I dynamics. Furthermore, there is a relatively sustained production of approximately 3.63 × 10 5 copies of genomic RNA from an estimated 5.37% of LTRs that persist in the intermediate state following activation (Fig. 2f). Since the presence of HIV-1 Tat protein is critical for activated transcription, we also investigated the dynamics of Tat expression as a function of time. As the results in Fig. 3g suggest, the production of Tat increases linearly over the modeled time frame measuring 5.19 × 10 5 proteins at 120 h indicating continual transcription activation resulting in the production of virus. To confirm this, we evaluated the production of genomic RNA (env) from a fully activated LTR, which we have modeled to elicit the production of fully infectious virions, over www.nature.com/scientificreports www.nature.com/scientificreports/ 120 h. The results in Fig. 3g confirm this showing an increase in production of genomic RNA. Following induction there is a steady increase in the number of genomic RNA transcripts with approximately 1.48 × 10 8 copies at 120 h. Altogether, our model effectively demonstrates the sustained production of small non-coding RNA, regardless of the activation state of the LTR, and the complex dynamics of genomic RNA production as a function of LTR activation.

LTR activation-dependent protein expression in T-cells.
Following transcription, the Gag polyprotein (Pr55) is expressed from unspliced genomic RNA. During maturation, Pr55 is cleaved by a virally encoded protease to produce four smaller viral proteins, including the HIV-1 capsid protein, p24, a marker of virus maturation. Along these lines, we examined the production of both Pr55 and p24 in T-cells following incubation in low serum media (0.1% FBS) for 36 h, subsequently cultured in 20% FBS media, and treated with an inducer for 24 h. Whole cell extracts were analyzed by Western blot using a p24 antibody. Not surprisingly, the results in Fig. 4a show a high level of p24 background expression in T-cells at 0 h (lane 1), which was expected due to higher levels of basal transcription in these cells. Furthermore, there is little to no detection of Pr55, suggesting efficient cleavage of the polyprotein to produce p24. Analysis of whole cell extracts at 24 h after introduction of an inducer and the addition of 20% serum media ( Fig. 4a; lane 2) indicate an increased production of Pr55 as compared to 0 h samples which subsequently results in the increased production of p24 mediated by viral protease cleavage of the polyprotein.
To determine parameter values for our mathematical model, densitometry counts were taken of the Pr55/ p24 Western blot (Fig. 4a) to estimate rates of protein production (α p2 ) and subsequent cleavage (α p3 ). As it can be assumed that the translation of Pr55 from genomic RNA is unaffected by the introduction of an inducer, α p2 has equivalent production rate values from both pools of genomic RNA, termed env I and env A (Fig. 1). To better represent the system, α p3 , the rate of Pr55 cleavage into p24, has been derived from densitometry counts corresponding to the increase in p24 production over 24 h, rather than the increase in p24 relative to Pr55 production. Additionally, although Fig. 4a demonstrates the presence of p24 at 0 h, we hypothesized that this is unlikely to be production of novel RNAs produced by transcriptional activation coming from LTR R . To test this, T-cells were cultured in low serum media (0.1% FBS) for 36 h followed by the addition of 20% FBS media. Cells were then lysed using a non-detergent lysis buffer and the lysates and conditioned supernatants were used to treat uninfected T-cells and myeloids +/− an entry enhancer. The results from the infectivity assay in Supp. Fig. 2a show no p24 in recipient cells when treated with lysates or supernatants from transcriptionally latent cells (0 hr; lane 3 and 7, respectively) in the absence an of entry enhancer. Interestingly, when an enhancer is added, p24 is present in recipient cells (lane 4 and 8). However, this p24 is unlikely to represent the production of virus as there is a distinct lack of gp120, the HIV-1 envelope protein. These results suggest that incubation in low serum media effectively induces a transcriptionally repressed state without the production of infectious virus and that the presence of p24 in recipient uninfected cells is likely due to the uptake of Gag-coding mRNA or genomic RNA from the cell lysates either as free RNA or RNA packaged into extracellular vesicles (EVs) 33,39,48 . Furthermore, these results confirm the production of fully infectious virus at 24 h from LTR A (Supp. Fig. 2a, lanes 5-6 and 9-10). These conclusions were confirmed using RT-qPCR (Supp. Fig. 2b-e). For clarity of the proposed model, 0 h p24 values have been set to zero.
The results in Fig. 4b show the modeled dynamics of Pr55 production over 120 h, with no Pr55 associated early time points as 100% of the LTRs are in the repressed state. This is followed by a slow increase in Pr55 production over the modeled time frame resulting in approximately 1.17 × 10 7 Pr55 at 120 h. Similarly, the expression of p24 shows little fluctuation over the modeled time frame (Fig. 4c). These results suggest there is a slow, but steady increase beginning at 0 h which is sustained throughout 120 h, potentially due to accumulation of Pr55 as the LTR switches from an intermediate state to one of full activation, resulting in 5.67 × 10 5 p24 at 120 h. Taken together, these results show an increase in the production of both viral proteins, Pr55 and p24, as a result of viral activation.

Modeling dynamics of LTR activation in macrophages.
To determine transcriptional differences between HIV-1 infected T-cells and myeloid cells, the same experiment was performed as described previously using U1 whole cell extracts for in vitro kinase labeling with [γ-32 P]-ATP with Histone H 1 as a substrate. Results in Fig. 5a show that, similar to T-cells, the presence of low serum media elicited low kinase activity levels at 0 h (lanes 1 and 4). When cells were introduced to 20% FBS media there was an observed increase in kinase levels (6 h, lane 2) followed by a reduction in activity at 24 h, although kinase activity levels did not return to basal levels (lane 3) as seen in T-cells (Fig. 2a, lane 3), suggesting sustained transcriptional activation in myeloid cells from the intermediate state of the LTR as compared to T-cells. Upon induction, myeloid cells followed a similar pattern with an increase in kinase activity at 6 h, however, there was a faster return to baseline levels at 24 h in comparison to T-cells (Fig. 5a, lane 6). These results indicate distinct transcriptional timing in myeloid cells as compared to T-cells with myeloid cells having a more sustained level of LTR I transcription.
Similar to the T-cell experiments, CDK9/Cyclin T 1 complex was immunoprecipitated from infected U1 whole cell extracts from 0 and 6 hrs (+/− inducer) and used for an in vitro kinase assay. The results in Fig. 5b, show low background levels of CDK9 kinase activity at 0 h, similar to J1.1 (Fig. 2b). Following introduction of 20% FBS media, there was a 55.5% increase in CDK9 activity at 6 h (lane 2) with a similar increase activation (lane 3; 54%) in the presence of an inducer. Taken together, these results demonstrate that basal p-TEFb transcription in HIV-1 infected myeloid cells can be increased upon the introduction of serum and that the activation from basal transcription in myeloid cells may be distinctly different from T-cells, potentially due to differences in transcriptional machinery.
To mathematically model the changes in the HIV-1 LTR dynamics in myeloid cells, parameter values were constructed as previously described. The results in Fig. 5c indicate the proportion of myeloid cell LTRs in a repressed state over the modeled time frame (120 h). These data show that there was in a sharp decline in the Scientific RepoRtS | (2020) 10:2227 | https://doi.org/10.1038/s41598-020-59008-0 www.nature.com/scientificreports www.nature.com/scientificreports/ relative amount of LTRs in a repressed state following activation of the HIV-1 LTR, resulting in the repression of only 2.77% of LTRs at 120 h. Figure 5d shows the relative proportion of intermediate state LTRs in U1 cells which result in the production of HIV-1 RNAs via basal transcription. Initially following activation with 20% serum, there was a drastic increase in the proportion of LTR I which resulted in a peak at 12.43 h post induction yielding approximately 37.80% of HIV-1 LTRs in an intermediate state. The observed increase was followed by a subsequent plateau (20.81%) which was sustained throughout the modeled time frame. Not surprisingly, there was an associated increase in the activated LTR state (LTR A ; Fig. 5e) following induction which resulted in approximately 76.42% of HIV-1 LTRs producing full length, genomic HIV-1 RNA and, in turn, replication-competent virions at 120 h. Altogether, the changes in LTR activation states post-induction in myeloid cells is shown in Fig. 5f. Although the overall patterns of change in LTR states was similar to that of T-cells, these data indicate LTR dynamics unique from that of HIV-1 infected T-cells in terms of timing, relative proportions of LTR states during extended time frames (120 h), and notably less variability (i.e. noise) as indicated by the confidence intervals.

Production of viral RNA in macrophages.
To compare production rates and copy numbers of short, non-coding and full-length HIV-1 viral RNA transcripts, infected U1 cells were cultured as previously described for T-cells. RT-qPCR of three biological replicates analyzed in triplicate was used to examine the production of two HIV-1 viral transcripts, TAR and env (i.e. genomic) RNA. For RNA production from LTR I , samples taken at 0 and 24 h in the absence of an inducer were analyzed and results in Fig. 6a indicate a high level of background transcription of TAR and genomic RNA associated with 0 h (4.70 × 10 6 and 2.96 × 10 5 , respectively). This was likely due to low levels of basal transcription during a latent state, a finding that is consistent with several in vivo studies 6,32,44 . Upon introduction of 20% FBS media, there was little increase in the transcription of TAR (6.27 × 10 6 copies) and genomic (4.91 × 10 5 copies) RNAs at 24 h.
RT-qPCR was also performed in technical triplicate to assess for viral RNA production in three biological replicates at 0, 24, 48, 72, 96, 120 h to determine production at LTR A . Data in Fig. 6b show a high level of background www.nature.com/scientificreports www.nature.com/scientificreports/ transcription of both RNAs at 0 h (3.13 × 10 6 and 1.49 × 10 5 , respectively). This was likely due to high levels of RNA prior to transcriptional suppression by low serum media. Upon introduction of 20% FBS media and PMA/ PHA, there was an increase in both TAR and env RNA in all measured time points relative to 0 h, with a peak in transcripts at the 72 h time point. Collectively, these data show an increase in the transcription of both measured RNAs following introduction of 20% serum media and that myeloid cell transcription is better activated by the introduction of an inducer (PMA/PHA).
To mathematically model our experimental results, we constructed RNA production rates as previously described. Similar to T-cells, we reasoned that the persistent, non-processive transcription observed in several studies involving in vitro and in vivo analysis of RNA production during antiretroviral therapy suggests that there is TAR RNA production during the repressed state of the HIV-1 LTR 6,44,47,48 . As such, we have estimated the production of short, non-coding TAR RNA in infected myeloid cells to happen at a rate of 2.9 × 10 4 copies/mL/h (α m1R ). As previously described, rates for TAR (α m1I ) and env RNA production (α m2I ) from LTR I were determined using experimental data from Fig. 6a. LTR A RNA production rates (α m1A and α m2A ) as well as RNA degradation rates (γ m1 and γ m2 ) were determined using parameter fitting (Fig. 6c,d). Densitometry counts of Tat Western blots (Supp. Fig. 1) were used to calculate a Tat production value (α p1 ).
The results in Fig. 6e demonstrate that the production of TAR RNA occurred in HIV-1 infected myeloid cells regardless of the state of the LTR, as evident by the background level of TAR RNA (3.13 × 10 6 copies). Furthermore, the production of short, non-coding RNA continued to increase throughout the measured time frame resulting in approximately 1.09 × 10 4 copies at 120 h. The production of env (genomic) RNA from the Whole cell extracts were then analyzed using an in vitro total kinase assay (a) or a CDK9 IP kinase assay (b) to evaluate for changes in the HIV-1 LTR. Biochemical assays were analyzed by densitometry, counts were normalized and used to construct parameters for mathematical modeling to assess relative proportions of the HIV-1 LTR in the various states; repressed (c), intermediate (d), and activated (e) over 120 h. The black line demonstrates the solved value of the original parameter set, grey lines are realizations due to the sampling of parameters using a Latin hypercube sampling method. The dashed green, red and blue lines represent 80%, 90% and 95% confidence intervals, respectively. (f) All three LTR states; repressed (LTR R , red), intermediate (LTR I , blue), activated (LTR A , green) are superimposed to indicate changes relative to each state at a given time. Selected lanes from the same blot with identical exposure settings are presented in panel a.

Scientific RepoRtS |
(2020) 10:2227 | https://doi.org/10.1038/s41598-020-59008-0 www.nature.com/scientificreports www.nature.com/scientificreports/  Fig. 6f. These data indicate a slight increase in the production of full-length genomic RNA from the intermediate state of the LTR over 120 h which yielded an estimated 3.03 × 10 5 copies from approximately 20.81% of the HIV-1 LTRs which maintained an intermediate state following activation. Similar to T-cells, we have also modeled the production of HIV-1 Tat protein. Figure 6g shows a modest increase in the production of Tat resulting in approximately 1.09 × 10 4 proteins at 120 h, indicating sustained activated transcription over the modeled time frame and resulting in the increased production of infectious virions. In agreement, Fig. 6h shows an increased production of env from an activated LTR state (LTR A ). Initially (0 h), 100% of the HIV-1 LTRs were in a repressed state due to the presence of low serum media. Upon addition of 20% serum media there was a continued increase in the production of env which resulted in approximately 2.14 × 10 6 copies at 120 h. Taken together these results indicate the continuous production of TAR RNA from HIV-1 infected myeloid cells and demonstrates the production of viral RNA as a function of the HIV-1 LTR activation state. Furthermore, these data suggest the magnitude of production of viral products may be inherently different in myeloid lineage cells as compared to T-cells.

Viral protein production in HIV-1 infected macrophages. To complement the T-cell experiments,
we analyzed the production of Pr55 (Gag polyprotein) and p24 in myeloid cells post-incubation with 0.1% FBS media for 36 h and subsequent treatment with an inducer for 24 h. Following treatment, cells were lysed and analyzed for the presence of these viral proteins using a p24 antibody. The results in Fig. 7a show very low overall levels of Pr55 and p24 background expression in U1 cells at 0 h (lane 1), which is in line with infectivity assay results in Supp. Fig. 2. Interestingly, these levels were markedly reduced as compared to T-cell whole cell extracts at 0 h ( Fig. 4a; lane 1). Following treatment with an inducer and the addition of 20% FBS media (Fig. 7a; lane 2) there was an increased production of Pr55 and p24, albeit expression was lower than that of J1.1 (Fig. 4a; lane 2).
To model these biochemical results mathematically, parameter values were constructed as previously described. Figure 7b,c indicate the changes in viral protein production over 120 h including the presence of no Pr55 (Fig. 7b) and p24 (Fig. 7c) at 0 h due to 100% of the HIV-1 LTRs in the repressed state (LTR R ). The data in Fig. 7b show a slow, yet steady, increase in the production of Pr55 following activation of the LTR by 20% serum and introduction of an inducer, which resulted in approximately 2.17 × 10 5 Pr55 at 120 h. The same trend was observed for p24 levels which yielded an estimated 5.84 × 10 3 p24 proteins (Fig. 7c). Not surprisingly, the data suggests overall lower levels of Pr55 and p24 in myeloid cells as compared to T-cells at extended time points. Interestingly, the 2-log difference between cell types is maintained for both Pr55 and p24 suggesting similar rates of Pr55 polyprotein cleavage in both cell types.

Comparing transcriptional differences in HIV-1 infected T-cells and macrophages.
To address differences in LTR state dynamics amongst cell types, solved parameter values for myeloid and T-cells were superimposed to determine timing and magnitude differences between cell types. The data in Fig. 8a demonstrates that myeloid LTRs exhibit faster changes in LTR states as indicated by rapid exit of the LTR from the repressed state, complemented by quick entry into an activated LTR as compared to T-cells. These cell-type differences are supported by a 2018 study that utilized a multiply spliced transcript reporter assay to measure the generation of RNA over time and found that primary myeloid cell multiply spliced transcripts, such as those which would be generated by LTR I , peaked more quickly than T-cells 50 . This trend was validated in primary latent peripheral blood mononuclear cells (PBMCs), which were separated into T-cell and myeloid populations and measured for rates of LTR state activation by the same methods as previously described for our cell lines (Supp. Fig. 3). Furthermore, these data show a more heterogeneous population of LTR states in myeloids given that, at extended time frames  Fig. 8b), likely due to the larger proportion of T-cell LTRs in the activated state. Tat production, produced from LTR I in our model, closely follows a similar trend (Fig. 8c), due to a lower calculated production rate for myeloids despite a larger proportion of myeloid LTRs in an intermediate state at 120 h. The production of full-length genomic RNA from both LTR I (green) and LTR A (blue) in T-cells and myeloid cells is shown in Fig. 8d. As expected, the production of env I (green lines) was higher than env A (blue lines) in both T-cells and myeloids at early time points. Interestingly, the transition from a predominately env I state to a predominately env A state occurred earlier in T-cells (approximately 28.5 h) than myeloid cells (approximately 50 h) despite the more rapid change in LTR state exhibited by myeloid cells. This is likely due to the high rate of RNA production in T-cells as compared to myeloids. Finally, the production of HIV-1 viral proteins, Pr55 and p24, was increased in T-cells as compared to myeloid cells (Fig. 8e). Collectively, these results suggest that T-cells produce viral products at higher levels compared to their myeloid counterparts although the myeloids show a more rapid change in the LTR state.
Modeling the use of F07#13. The dependence of activated viral transcription on the HIV-1 Tat protein makes it an attractive target for therapeutics aimed at inhibiting transcription of the HIV-1 virus. As a result, Tat peptide mimetics, specifically F07#13, have been developed and shown to inhibit transcription both in vitro and in vivo 15,51 . To further test our model, we tested transcriptional changes associated with the introduction of the transcriptional inhibitor F07#13. F07#13 treatment was incorporated into the governing equations of the model within the terms k A (Tat) and k I (Fig. 1) . Data in Fig. 9a  www.nature.com/scientificreports www.nature.com/scientificreports/ reduction in p-Histone H1 (i.e. maximal drug effect) in HIV-1 infected myeloids at 6 h post-induction (U1; lanes 2 and 5), suggesting differences in transcription between cell types and variation in drug efficacy. The simulation data in Fig. 9b,c show the mathematically derived prediction of the relative proportions of LTR states ± F07#13 in T-cells and myeloids, respectively. As described in Figs. 2f and 5f, in the absence of F07#13, there was a rapid decrease in the number of LTR R (red line) concomitant with a sharp increase in the proportion of LTR I (blue line) until ~21.31 h and 12.43 h in T-cells and myeloids, respectively. Furthermore, there was an associated increase in activated LTRs (green line) following induction. The graphs show a decrease in the proportional of LTRs in the active state in response to F07#13 treatment (dotted green line). Interestingly, there was a slight increase (9.35%) at 120 h in the relative amount of LTRs in the basal transcriptional state (LTR I ; dotted blue line) observed in T-cells (Fig. 9b) which is evident at approximately 9 h and was sustained through the modeled time frame. In myeloid cells there was a larger increase (35.66%) at 120 h in the relative proportion of LTR I , evident at approximately 6 h, which was similarly sustained throughout the 120 h time frame. Finally, there was a slight increase in the number of LTRs in latency (LTR R ; dotted red line) in both cell lines, suggesting preferential basal transcription during F07#13 treatment and thereby the production of non-coding RNAs such as TAR. In agreement with our biochemical assays, graphs generated from our model demonstrate the F07#13-mediated LTR effect in T-cells peaks at ~24 h (Fig. 9a) whereas the same effect was observed earlier in myeloid cells (6 h; Fig. 9a). To illustrate the potential utility of the model to explore the effects of F07#13 on the production of viral proteins and, in turn, infectious virions, the production of Pr55 in the presence of F07#13 was modeled and compared to the production of the protein in the absence of the drug. These data show the post-treatment outcomes for T-cells (Fig. 9d) and myeloid cells (Fig. 9e) in terms of viral protein production. Interestingly, these results indicate that F07#13 is more effective in lowering virus production in myeloid cells as compared to T-cells despite a decrease in the relative proportion of LTR A in both cell types. The decrease efficacy in T-cells is likely due to compensation by LTR I in the production of env resulting in no change in the production of Pr55. These data are in line with our previously published results which found that RNAPII was still present on the HIV-1 promoter post-treatment Western blots were assessed via densitometry, counts were normalized and used for calculation of protein production parameter values which were utilized for mathematical modeling of viral production of Pr55 (b) and p24 (c). The black line demonstrates the solved value of the original parameter set, while the grey lines are all the realizations with respect to the sampling of parameters using a Latin hypercube sampling method. The dashed green, red and blue lines represent 80%, 90% and 95% confidence intervals, respectively. Selected lanes from the same blot with identical exposure settings are presented in panel a. (2020) 10:2227 | https://doi.org/10.1038/s41598-020-59008-0 www.nature.com/scientificreports www.nature.com/scientificreports/ with F07#13 using an in vitro transcription assay 31 . Conversely, the addition of F07#13 in myeloid also decreases the relative proportion of LTR A , however, unlike T-cells this proportionate drop outweighs the relative increase in the proportion of LTR I .

Discussion
Here, we have generated a mathematical model to describe three states of HIV-1 transcription. This model was created in such a way that it is open to broad applications. We utilized 100% HIV-1 infected cell lines of two different immune cell types to show differences in the production magnitude and induction speeds between infected T-cells and myeloids. Specifically, our model demonstrates that latent infected myeloids become activated sooner than T-cells under the influence of a transcriptional inducer; however, latent infected T-cells generate more viral products than myeloids, despite a slower activation rate (Fig. 10). Interestingly, the rates of LTR activation in primary cells followed a similar trend to that seen with cell lines, although the overall rates were lower (Supp. Fig. 3). This could potentially be due to the fact that the primary cells were infected in vitro and likely not 100% infected. The fact that myeloid and T-cells show variances in the speed of LTR state change and levels of viral components produced points to differences in transcriptional machinery between the cell types. This notion has been supported by previous findings that have shown many differences in transcriptional machinery components, particularly host transcription factors, between myeloid and T-cells, which has been previously described at length by Rohr et al. 52 . The lack of an FDA-approved transcription inhibitor for HIV-1 illustrates the need for further research into compounds with the potential to limit viral transcription. As a result, we have integrated potential therapeutic interventions, such as F07#13, into our model to explain their effects on HIV-1 activation and generation of viral products.
The observed high production of TAR RNA in both cell lines (10 6 to 10 8 ) at all time frames as compared to the relative levels of other viral proteins and RNAs, suggests that regardless of cell type, HIV-1 produces short, non-coding RNAs (Figs. 3e and 6e) from all states of the LTR. These findings are confirmed by several patient and in vitro studies which show the presence of TAR RNA within infected cells despite antiretroviral therapy, which is effective in reducing other viral products and suppressing the virus in the plasma 6,32,44,48 . Several studies have previously shown that TAR RNA is released from HIV-1 infected cells, including during latent infection, which can then induce numerous inflammatory responses in neighboring cells 33,35,48 .
Other research suggests that during activated transcription, which results in the production of genomic RNAs, the ratio of TAR:Tat:Cyclin T 1 is 1:1:1 [53][54][55] . Furthermore, the presence of Cyclin T 1 greatly enhances the binding affinity of Tat to TAR, with a Hill coefficient of 2.7 and a dissociation constant (K D ) of 2.45 nM 56 . However, we propose that the production of TAR RNA as a short non-coding RNA within the cell could result in the binding of TAR RNA to Tat protein in the absence of Cyclin T 1 , potentially resulting in TAR RNA produced at high copy numbers to inhibit Tat-activated transcription. Along these lines, our experimental data show a relative abundance of TAR RNA in infected cells as compared to other viral products (i.e. Tat), which could potentially indicate that TAR RNA could be utilized to sequester Tat protein, thereby slowing activated viral transcription and promoting latency of the virus. Altogether, this finding may potentially point to a novel mechanism of viral maintenance of latency.
It has been previously shown that TAR RNA is released in EVs such as exosomes both in vitro and in vivo even under antiretroviral regimens 33,35,44,48 . Therefore, the γ m1 and γ m2 parameters include exit of viral products from the cell through these secretory mechanisms. Additionally, the fact that there is a high levels of TAR RNA in T-cell EVs may suggest enhanced efficiency of T-cells in the packaging of EV products. The data in Supp. Fig. 5 support this hypothesis and demonstrate that the EV-associated viral RNA copy number (white bars; TAR, Supp. Fig. 5a; env, Supp. Fig. 5b) is higher in T-cells as compared to myeloids. Furthermore, there are innate differences in the distribution of viral RNA (TAR and env) within T-cells and myeloids. T-cells exhibit approximately equal copy numbers of both TAR (Supp. Fig. 5a) and env (Supp. Fig. 5b) within the nucleus (black bars) and cytoplasmic fractions (grey bars) while myeloids show a higher proportion of TAR (Supp. Fig. 5a) and env (Supp. Fig. 5b) in the cytoplasmic fraction in comparison to the nuclear fraction. These data may suggest differences in the nuclear export of viral RNA which could potentially lead to variations in protein expression as previously demonstrated by the model. Additionally, these data imply that a TAR-mediated Tat inhibition could potentially play a larger role in T-cell latency as compared to myeloids. These observations reiterate broad utility of the model shown An interesting finding within our biochemical assays showed that there were background levels of p24 consistently produced from latent T-cells suggesting persistent production of mature Gag. Previously, we have shown that exosomes contribute to the lack of transcriptional latency within HIV-1 infected cells 32 . Therefore, the presence of p24 despite a serum-starvation induced latent state further supports the notion that there may not be a true transcriptional latency during HIV-1 infection of T-cells (Fig. 4a). This is in contrast with our observation in myeloid cells (Fig. 7a) which showed minimal levels of p24 in serum-starved cells suggesting different mechanisms of maintaining viral latency between these cell types. This could pose a potential hurdle in designing targeted antiretrovirals and the development of future transcription inhibitors against HIV-1.
The mathematical model presented in this paper is novel in that it incorporates new biological mechanisms in a way that preserves important biologically realistic features. Work is currently under way to study analytical properties of this system, and several important facts have been already established. For example, starting with a positive amount of LTR R , positivity of all variables is maintained during evolution. Furthermore, all variables are bounded, which is consistent with the fact that no realistic environment can support indefinite exponential growth. During the construction of the model, some simplifying assumptions were made regarding the switching function, A (Tat) and 2 (Tat), and the production/decay parameter. This is reasonable because, for a short time frame (5 d), step-wise constant rate approximation is sufficient for capturing the complexity of the system behavior. However, we speculate that nonlinear forms of the switching functions and more biologically realistic forms of the rates of change would be more suitable to study the long-term behavior of the system. While simplistic, the linear system allows for a more thorough understanding of the transient modes, for instance, it may be helpful in clarifying how different rates contribute to the dynamics of each variable at the beginning of transcription. In this manuscript, parameter estimations were carried out for a relatively small subset of parameters that cannot be directly inferred from the data. In future work, a more thorough sensitivity analysis can be performed to establish the role of each of the parameters.
In the future, additional cell types, including primary cells and other immune cell types, or various drug treatments designed to act at multiple stages of latent virus reactivation may be easily integrated into the model we present here to show their downstream effects and optimal time frame for usage. To demonstrate this utility, we have integrated the transcription inhibitor F07#13 into the model to show elicited changes in LTR activation post treatment, which varied between cell types. Specifically, F07#13-repressed transcription was more quickly achieved in myeloid cells in comparison to T-cells (Fig. 9b,c). Furthermore, the differences in the relative proportions of LTRs during various time points may dictate the efficacy of the drug in terms of virus production (Fig. 9d,e). One of the potential applications of the proposed model is the study of the effectiveness of drugs in silico. While this type of work can be carried out numerically on a case-by-case basis, we suggest that valuable insight may be derived by first studying the model analytically and establishing conditions when a certain drug is effective, especially when it is used in combination with other drugs. Other applications of the dynamics of this system or its nonlinear version may potentially be extended to other retroviruses including endogenous ones that are known to activate under given stimuli. An example of this could be the activation of human endogenous retrovirus type-K (HERV-K) by HIV-1 Tat 57 . Moreover, the model could be used to more closely study the various states of the HIV-1 LTR (repressed, intermediate, and activated). For the repressed state of the LTR, these studies www.nature.com/scientificreports www.nature.com/scientificreports/ could include exploring the effects of various suppressive transcription factors, transcriptional gene silencing mechanisms, and suppressive chromatin complexes on relative proportions of the LTR transcriptional state and the production of viral RNAs and proteins. Inquiries into the intermediate state of the LTR such as RNAPII phosphorylation, basal transcription factors, chromatin looping, and enhancers which could also be incorporated into the model. Additionally, Tat-activated transcription studies such as RNAPII elongation through various genome regions and RNA processing (capping, splicing, 3′end processing) could also be explored. For example, Yukl et al. has designed a set of primer and probe sequences which could be easily incorporated to increase the resolution of the model presented here by analyzing several different RNA transcripts. This approach could aid in fine-tuning the model as it includes the use of droplet digital PCR to measure transcripts associated with various blocks in transcription which would allow for more precise inclusion of drug mechanisms 58 . Finally, the presented model could be used to larger data sets regarding integration sites of the HIV-1 genome or other viral genomes and the downstream effects in terms of transcriptional regulation of viral products using chemical or other biological compounds that serve as therapeutics.

Materials and Methods
Cell culture and reagents. U1 (HIV-1 infected promonocytic cell with 2 copies of virus; one wild-type and one mutant) and J1 (HIV-1 infected T-cell with one copy of wild-type virus) cells were cultured and maintained in RPMI 1640 supplemented with 10% fetal bovine serum (FBS), 1% penicillin/streptomycin, and 1% L-glutamine (Quality Biological) and incubated in 5% CO 2 at 37 °C. U1 cells were activated with phorbol 12-myristate 13-acetate (PMA; 1 µM) and phytohemagglutinin (PHA; 10 µg/µL) and then treated with F07#13, an allosteric Tat/ CDK9 inhibitor (1 µM), or irradiation (5 Gy). Treated cells were incubated in 5% CO 2 at 37 °C. Both cell lines (U1 and J1) were provided by the AIDS Reagent Program (National Institutes of Health). Both U1 and J1 cells were serum starved for 3 days, washed, counted, and plated in 20% FBS media +/− inducer. Samples were collected at 0, 6, 12, 24, 48, 72, 96, and 120 h post-induction. The 0 h sample represents samples that were collected 10 minutes post-treatment due to the size of experimental design and handling procedures. Therefore, these samples would have shown some level of transcription within the first 10 min, since the rate of RNA Pol II transcription is fairly rapid and high (the RNAPII complex moves through very long introns and also through regions dense with alternating exons and introns at an average rate of ∼3 kb per min) 59 . Peripheral blood mononuclear cells (PBMCs) were treated with PHA/PMA and IL-2 treatment. After one week, they were infected with HIV-1 (89.6; multiplicity of infection = 3.0). Adherent macrophage cells were separated from T-cells using centrifugation, and 3 days post-infection, T-cells were treated with cART and IL-7, which allowed cells to enter latency. Adherent macrophages were also treated with cART to obtain latent cells. Prior to activation with 20% FBS media and PMA/PHA, cells were kept in low serum media to suppress basal transcription.
Kinase assay. In vitro kinase assays were carried out using whole cell extract (prepared as described below). www.nature.com/scientificreports www.nature.com/scientificreports/ the BioRad CFX96 Real Time System. Quantification of the samples was determined based on the cycle threshold (Ct) value relative to the standard curve.

Preparation of whole cell extracts and western blot analysis. Infected cell pellets were harvested
and washed with 1X PBS without calcium and magnesium. The resulting pellet was resuspended in lysis buffer [50 mM Tris-HCl (pH 7.5), 120 mM NaCl, 5 mM EDTA, 0.5% Nonidet P-40, 50 mM NaF, 0.2 mM Na 3 VO 4 , 1 mM DTT, and 1 complete protease inhibitor mixture table/50 mL (Roche Applied Science, Mannheim, Germany)]. The mixture was incubated on ice for 20 min with vortexing every 5 min followed by separation via centrifugation at 10,000 × g at 4 °C for 10 min. Lysate total protein concentration was assessed using Bradford protein assay (Bio-Rad).
Isolation of nuclear and cytoplasmic compartments. Cytoplasmic and nuclear compartments were separated and extracted using the NE-PER Nuclear and Cytoplasmic Extraction Reagent kit (Thermo Fisher Scientific) as per the manufacturer's instructions. Briefly, log-phase J1.1 and U1 cells were collected and resulting culture supernatant was used for exosome isolation as described below. Cell pellets were washed in 1X PBS and incubated with 100 µL of cold CER I followed by vortexing (10 seconds) and subsequent incubation on ice for 10 min. Following incubation, 5.5 µL ice cold CER II was added to both samples and vortexed for 10 seconds followed by a 1 min incubation on ice. Samples were then vortexed, centrifuged for 5 min at 16,000 × g to separate the cytoplasmic compartment which was transferred to a new, pre-chilled tube. Pelleted nuclei were treated with 50 µL of ice cold NER, vortexed for 20 seconds and incubated on ice for 45 min with additional vortexing (every 10 min). Following incubation, samples were centrifuged at 16,000 × g for 10 min and resulting nuclear extract was transferred to a new prechilled tube. Both cytoplasmic and nuclear fractions were then subjected to subsequent RNA isolation, reverse transcription, and qPCR as described above.
Isolation of exosomes. J1 and U1 culture (five days post last feed; late log phase of growth) supernatant was harvested and enriched for exosomes using Nanotrap particles (NT80/82) as previously described [32][33][34][35]48,60 . Briefly, 30 µL of a 30% slurry of NT80/82 (1:1) was added to 1 mL of cell culture supernatant and incubated with rotation overnight at 4 °C. The next day NTs were pelleted and washed with 1 mL of sterile 1X PBS and subjected to RNA isolation, reverse transcription, and RT-qPCR as described above.
Densitometry. Raw densitometry counts were obtained using ImageJ software. Exposures were matched between membranes according to positive control signals. Densitometry data was normalized by a two-step process. First, the background measurements for each membrane were subtracted from each band of interest. Second, each protein band's raw densitometry count was converted to a percentage by dividing it by the background measurement for each membrane and multiplying by 100%. This normalization accounted for differences in the backgrounds and exposures between membranes, while also converting densitometry counts to an appropriate unit for use in our model. Statistical analysis. Standard Deviations (S.D.) were calculated from technical triplicates for every quantitative experiment using Microsoft Excel. To determine the sensitivity of the model is to its parameters, we employed the Latin Hypercube Sampling method. First, we assumed each parameter is normally distributed, whose mean is the estimated value from the data (either experimentally and numerically derived) and standard deviation is taken to be 1% of the parameter value. Using this method, we selected 200 parameter sets and generated 200 realization curves for the model, which are plotted along with the solution curve. The 80%, 90% and 95% confidence intervals were constructed to contain 80%, 90% and 95% of the realization curves, respectively.
Model. The mathematical model presented in this manuscript aims to describe key features associated with basal and activated transcription of the HIV-1 genome. In order to produce a model which closely represents the series of events which occur in vivo, the associated model parameters have been determined using in vitro experiments in HIV-1 infected T-cell (J1.1) and macrophage (U1) cell lines, with the exception of transcription rates from the activated promoter state (α m1A and α m2A ), which we obtained by numerical fitting of our model predictions to the experimental measurements of TAR and env. Our model, as shown in Fig. 1, categorizes the long-terminal repeat (LTR) into three states; a repressed state (denoted by LTR R ), an intermediate state (denoted by LTR I ), and an activated state (denoted by LTR A ), similar to that used in a previous report 27 . Furthermore, the relative proportions of viral RNAs, including short, non-coding RNA (TAR) and long genomic RNA (env), and resulting functional proteins have been attributed to each LTR state.
Activation of the HIV-1 LTR is dependent on Tat association with CDK9 and Cyclin T 1 . The phosphorylation of Histone H1 is directly correlated with exposure of the promoter and transcriptional activation. As such, a Histone H1 kinase assay (see Materials and Methods section) was used to quantitate activation of the LTR at the three transcriptional states (repressed, intermediate, and activated). The repressed and intermediate HIV-1 promoters are considered to be in the 'OFF' state while the activated HIV-1 promoter is considered to be in the 'ON' state. Transitions between these states occur with rates k ON Tat-dependent and has be described as a step-wise function. In this work we have used constant approximation for this parameter due to extremely rapid transition to the active state. There have been several reports of sustained intracellular levels of TAR RNA during antiretroviral treatment which can range from 1 × 10 1 to 1 × 10 5 RNA copies in vivo 6,44 . Therefore, assuming a 1:1 binding ratio of TAR:Tat, we conservatively estimated the critical Tat value to be equal to 3001 to allow a large enough pool of Tat protein to overcome potential sequestration by TAR RNA 33,35,48 , and thereby allow for efficient Tat-activated transcription. As a result, we propose that these LTR states and associated rates can be expressed as We have interpreted the three LTR states as proportions of the total LTR count so that LTR R + LTR I + LTR A = 100%. The rates k ON , k OFF , k A (Tat), and k I have units of inverse time with the change in each LTR state being a (unitless) proportion of the total LTR count.
A key new feature of our model is the inclusion of three different transcriptional states as described by in vitro quantitation of two different HIV-1 RNAs -long, genomic RNA transcripts which we have denoted as env I and env A (for transcripts that associate with the corresponding promoter states LTR I and LTR A , respectively)-and one short, non-coding RNA sequence which is denoted by TAR. HIV-1 TAR RNA is produced from a virally infected cell regardless of latency or suppressive cART therapy 33,35 . TAR generation occurs from each of the three LTR states (repressed, intermediate and activated) at specified rates. This short RNA sequence is transcribed by all three promoter states with associated rates as indicated by the following biochemical reactions: The production of Tat (i.e. translation of mRNA into a functional protein) occurs in our model through the intermediate, env I , long RNA. The production of Tat further enhances activation of the intermediate LTR state, giving a possible gene-activated feedback mechanism, as well as enhances transcription of activated long RNA transcripts env A , giving a transcription-activated feedback mechanism 61 . The expression of transcripts generates more virions that go on to infect more cells and continue the HIV-1 gene regulation. The Tat protein produced from the long, full-length transcript env I from the intermediate state of the LTR is expressed as Both long RNA transcripts, env I and env A , lead to the production of the HIV-1 Gag polyprotein, Pr55, which in turn produces p24. The levels HIV-1 Pr55 and p24 were visualized using Western blot analysis and quantitated using densitometry. Densitometry counts of Pr55 were then utilized to find the rate of Pr55 production (α p 2 ) from both env I and env A . The production of the HIV-1 capsid protein p24 which is cleaved from the gag polyprotein Pr55 at a calculated rate α p 3 . Finally, we expressed degradation and exit of viral products into the extracellular space through the parameters γ p 1 . (degradation of Tat), γ m 1 (degradation/exit rate of TAR), γ m 2 (degradation/exit rate of env I and env A ), and γ p 2 (degradation/exit rate of P 24 ). As we are primarily interested in the dynamics of these quantities over a short time frame (e.g. several days) simulations were run with the degradation/exit rate rates set to zero for Tat and p24. (2020) 10:2227 | https://doi.org/10.1038/s41598-020-59008-0 www.nature.com/scientificreports www.nature.com/scientificreports/ The first three of these differential equations (Eqs. 11-13) describe the evolution and transitions that occur between the three LTR promoter states. The dynamics of these three quantities depend on the dynamics of Tat through the activation rate k A (Tat). However, if the coefficient k A (Tat) is constant, then the LTR states are mathematically decoupled from other dynamic variables in our model. Our initial computations have used k A (Tat) as a constant although measured Tat-dependent rates can be easily incorporated. The three states LTR R , LTR I and LTR A are defined as proportions of the total LTR count and thus we have LTR R + LTR I + LTR A = 100%. That is, the sum of the three LTR state variables, which represent proportion of the total number of LTR sequences, is preserved. The fourth equation (Eq. 14) gives the rate of change of Tat which is produced by the long transcript env I and degrades at rate γ p 1 . The fifth equation (Eq. 15) characterizes the dynamics of the transcript TAR, which is generated at various rates by all three promoter states, LTR R , LTR I and LTR A . The sixth and seventh equations (Eqs. [16][17] describe the dynamics of intermediate and activated transcripts env I and env A which are generated by their counterpart promoter states LTR I and LTR A and are translated to Tat and Pr55 as well degraded. Finally, the eighth and ninth equations give the rates of change of Pr55 and p24. In addition to describing the dynamics of the transcription process, we have incorporated into the system of differential equations a means to assess the effectiveness of F07#13, Tat peptide mimetic. We have modeled F07#13's interaction with the HIV-1 LTR as an increase in the transition from LTR A to LTR I and simultaneous decrease in the transition from LTR I to LTR A . Thus, taking the effect of F07#13 into account, the rate of transition from LTR I to LTR A becomes k (Tat)LTR w A I 1 5 and the rate of transition from LTR A to LTR I becomes w 1 k 1 LTR A . Here, w 1 and w 5 are constants representing the effect of F07#13 which have value greater or equal to 1 that have been obtained experimentally. These constants have been incorporated into the previously described governing equations as follows: Parameter fitting. To account for the uncertainty associated with measurements of TAR and env in T-cells and macrophages, a numerical means was used to estimate four parameters: m1A , m2 , m1 and m2 . The model was fit simultaneously to the experimental values of TAR and total env (e.g. time series data of TAR and env performed in three biological replicate and analyzed in technical triplicate by RT-qPCR) independently for T-cells and Macrophages. In order to initialize the parameter estimation process, experimentally derived rates for the production of TAR and env were used as initial approximations of the parameter values. The initial approximations of the degradation/exit rates for TAR and env are taken as the linear rates of change of extracellular TAR and env RNA. The objective function to be minimized is the standard objective in least squared fitting. All time series data was simultaneously fit for each triplicates and replicates of both env and TAR. The function fmincon -a commercial product of MATLAB -was used to estimate these parameters within reasonable ranges. The ranges for the estimation extend from 0.1 to 10 times the initial approximations. Through testing, changing this range does not affect the overall results. The estimated values are given in Table 1. Impact statement. We propose a novel 3-state LTR model (LTR R , LTR I , LTR A ) incorporating long and short RNAs that offers insight into HIV-1 transcriptional dynamics in latently infected macrophages and T-cells and provides a framework for future mathematical analysis of HIV-1 transcriptional interactions and drug interventions. Specifically, presence of short non-coding RNAs such as TAR at high levels such as that produced from an activated LTR state may explain Tat sequestration and transcription inactivation leading to viral latency.

Data availability
The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.