1. Introduction
Cognitive impairment in HIV and dysregulation of dopaminergic systems are well-documented by us and others[
1,
2,
3,
4,
5,
6,
7,
8,
9,
10,
11,
12,
13,
14]. These clinical and physiological disturbances may result from substance use disorders that are frequent among people living with HIV (PWH)[
15,
16]. Psychostimulants such as methamphetamine (Meth) and cocaine are potent inducers of dopamine (DA) in the central nervous system (CNS). HIV enters the CNS during acute infection and remains latent or replicates at low levels in microglia reservoirs. While DA is fundamentally involved in mechanisms of reward, cells targeted by HIV in the CNS and elsewhere express all five DA receptors[
1,
11], with phenotypic characteristics that can be deeply affected by the DA-rich brain environment of substance users[
1,
2,
5,
7,
13,
17,
18,
19,
20,
21,
22,
23]. Indeed, many reports indicate that psychostimulants such as methamphetamine increase cellular susceptibility to HIV entry and spread[
1,
2,
11,
24], which ultimately may be critical for increased viral load in the brain of users as observed in animal models and human cohorts[
25,
26]. A question that remains, however, is the mechanistic basis for more frequently detectable viral load in brain and cerebrospinal fluid (CSF) from substance users and, in particular, what the role of DA in the increased susceptibility is. Currently, an outcome of effective antiretrovirals is the suppression of HIV in peripheral sites and fluids. However, latently infected long-lived innate immune cells with integrated viral sequences can either reactivate and reseed the infection or maintain a low-level production of HIV proteins with active roles in host cell chronic pathogenesis, even in PWH with viral suppression on ART [
27,
28]. We have previously shown that DA is a factor involved in cell-to-cell HIV spread[
1], but most importantly, in the latency reversal that is detectable in the CSF of methamphetamine users[
2]. Here, we report experiments in a promonocytic cell line latency model [
29] conducted to understand how elevated levels of DA found in the brains of methamphetamine users[
30] may trigger reactivation of the latent virus from the reservoirs of myeloid origin. The choice of U1 cells as a model was based on the consistency and predictability of the latent state without variables linked to technical factors with the goal of developing a strategy to facilitate the identification of virally infected/latent cells in situations of lack of reporters, or within the heterogeneity of myeloid target cells in vivo, preventing sorting or recognition.
Our previous studies have indicated that DA reverses latency in a significant fraction of, but not all, latent promonocytes[
2], leading us here to probe the use of single cell transcriptome analysis to compare features triggered by DA in cells that show signs of reversal in comparison to cells that maintain latency. The single-cell strategy enabled the detection of differences in the transcriptional landscape of all nuanced subpopulations showing elevated transcription of HIV gag-pol-env-tat-nef or no detectable viral transcription in response to DA, as well as its signatures, to define the requirements of latency regulation and pathways that may be harnessed for understanding neuroimmune interactions in populations of drug users. The results presented here further confirm the role of DA in latency reversal but also suggest exclusive pathways and orchestrated transcriptional programs that can be explored to control and shift the HIV latent state to active phenotypes.
This study specifically aimed to develop strategies that can help elucidate the mechanisms by which elevated DA levels trigger the reactivation of HIV from latency in myeloid cells and distinguish cells that remain latent for contrasting the heterogeneity of transcriptional responses of DA-sensitive and DA-insensitive cellular populations. By employing single-cell RNA sequencing to an in vitro myeloid latency model, we seek to identify unique gene expression patterns that may distinguish these cellular responses, thereby contributing novel insights into the neuroimmune interactions applicable to brain-derived targets pertinent to HIV persistence in substance users.
2. Materials and Methods
Promocytic Cell Lines and Culture Conditions – Chronically infected HIV-1 promonocytic (U1, ARP-165) cell lines were originally derived by limiting the dilution cloning of U937 surviving an acute infection with HIV-1[
29]. The U1 cells (ARP-165) [
29] were obtained from the NIH HIV Reagent Program, Division of AIDS, NIAID, NIH, generously contributed by Dr. Thomas Folks. The parent U937 cells (CRL 1593.2) were obtained from the American Tissue Culture Collection (ATCC, Manassas, VI). Both U1 and U937 cells were cultured in RPMI 1640 containing 2.0 mM of L-glutamine and 10% heat-inactivated fetal bovine serum (Lonza Bioscience, Morrisville, NC, USA), and maintained in the log phase with >98% viability before plating at 10
6/mL in 12-well plates and stimulation. The cells were used with no more than 5 passages.
Cell culture treatments – DA concentrations were optimized before experimentation. Dopamine hydrochloride (H8502, Sigma Aldrich, St. Louis, MO, USA) was used at 1 and 10 μM, to mimic median range levels found in the brain of drug users[
30] and did not produce significant differences in a pilot. Single cell results that are shown were obtained with 10 μM DA. Recombinant human heterodimer MRP8/14 proteins (S100A8/S100A9, R&D Systems, Minneapolis, MN, USA) were used at 1 μM, as previously described[
2].
Single cell RNAseq samples preparation– Single cell RNAseq were performed using 10x Genomics Chromium Next GEM Single Cell 3ʹ Kit v3.1 (10X Genomics, Pleasanton, CA) as per manufacturer recommendations. Briefly, upon stimulation, the cells were collected by centrifugation (400 x g, 5 minutes, at room temperature), and washed with 1 ml PBS containing 2% of bovine serum albumin (BSA, Sigma-Aldrich). We loaded 10,000 fresh cells along with gel beads containing unique barcodes and reagents, on 10X Genomics Chip G, and using the 10x Chromium Controller that perform the gel-emulsion (GEM) single-cell droplets mRNA capture process. Within each GEM, cell lysis occurred, and mRNA was copied into cDNA by reverse transcription, tagged with a cell barcode and UMI (Unique Molecular Identifier). The GEMs were then broken and cDNA was cleaned up and amplified. Using Chromium single cell 3’ Reagent kits with v3.1 Chemistry Dual Index reagents, the cDNA was subjected to fragmentation and Illumina sequencing-ready library by the addition of sequencing adapters. The size of libraries as well as quantities were respectively measured by using Agilent 2200 TapeStation System (Agilent, Santa Clara, CA) and QubitTM 1x dsDNA HS Assay Kit (Invitrogen, Q33231) before testing. A Single Cell 3’ Gene Expression library was generated using a 25% of the total extracted cDNA. Samples in duplicates were sequenced on Illumina NovaSeq 6000 with an average depth of at least 20,000 100bp x 100bp PE reads per cell.
Single cell RNAseq data analysis
(i) Alignment and Quality filtering – Sequencing data processing were performed using Cloud-based Cell Ranger 7.2.0 platform and visualized in Loupe 7.0 browser (10X Genomics). Specifically, the scRNA data was demultiplexed and aligned to the customized hg38 Homo sapiens reference genome (NCBI RefSeq assembly GCF_000001405.40) which was combined with a chromosome representing the HIV-1 genome (AF033819). The generated sequencing libraries were then aggregated using the function Cell Ranger Aggr (v3.1.0) and low quality single-cell transcriptomes were filtered based on UMI count [500 to 30,000], gene count [300 to 5,000], and mitochondrial percentage [less than 15%]. Table 1 shows Cell Ranger derived quality controls.
(ii) Data normalization, batch correction, and clustering analysis. Following data normalization in Seurat v4 (
https://satijalab.org), batch correction was performed using Harmony. Using the 2000 most variable genes (default parameters), and 50 most significant principal components, we performed linear dimensional reduction and built neighborhood graph using the t-SNE low-dimension visualization tools in Seurat. Cluster-specific pathways and biological processes were analyzed using iPathwayGuide (AdvaitaBio, Ann Arbor, MI) [
31], as well as the Database for Annotation, Visualization and Integrated Discovery (DAVID)[
32], and visualized using GeneMania[
33] in Cytoscape 3.10.2[
34,
35].
(iv) Differential gene expression (DGE) analysis was done with DESeq2 R Package. We used the Cell ranger integrated to calculate DGEs and generate volcano plots for visualization of differential gene expression upon different stimulations: non-treated (NT) vs. Dopamine (DA), within and between U937 and U1 cells (Table 1). Statistical significance was tested using the Wilcoxon matched pair signed rank test. Cell Ranger data was converted to Seurat for aggregation and cluster analysis.
(v) All raw and processed data were deposited in GEO with the assessment number GSE278043.
Table 1 - Cell Ranger-derived quality assessment and count summary.
Systems Biology and Visualizations
Each cluster was treated as an individual experiment for input in .txt format, used to identify significantly impacted signaling pathways, Gene Ontology terms, disease processes, predicted upstream regulators, and putative mechanisms based on significant gene expression signatures. Comparisons between clusters of interest were performed using meta-report features in iPathwayGuide. Gene network analysis was performed in Cytoscape v3.10.2 [
35], with GeneMania plug-in[
36,
37], and Homo sapiens sources from BioGRID_ORGANISM[
38,
39,
40]. Permutation analysis in iPathwayGuide was used to rank gene lists by cluster matched to Gene Set Enrichment Analysis (GSEA) database[
41,
42] using keyword “HIV”.
3. Results
The transcriptional landscape of latently infected U1 cells and their parent uninfected U937 cell line was compared in the absence and presence of DA treatment. The integration of heterogeneous features in the dataset was visualized in a probabilistic framework using uniform manifold approximation and projection (UMAP) for dimentional reduction (
Figure 1). DA or vehicle exposure were carried out for 24 hrs, and the cells were prepared for single-cell RNAseq as described above. We identified a marked heterogeneity within the cell lines at any given condition and shifts due to treatment, detectable in unsupervised analysis both by dimensionality (
Figure 1A) and by subset clusters with common and exclusive profiles (
Figure 1B and 1C). The data revealed that the exposure to DA deeply modified gene expression and signatures, both in uninfected U937 and in latent U1 cells (
Figure 1A). Dimension reduction visualization of individual and merged conditions indicated cell clusters with exclusive transcriptional profiles resulting from HIV latent virus (clusters 2 and 10) or DA (clusters 1 and 8), as well as their interactions (clusters 3 and 6), (
Figure 1B and 1C), and clusters with overlapping dimensions across conditions (See Supplementary Table 1 for detailed assignment of clusters to experimental conditions). The observation of individual clusters with overlapped single cell events associated with transcription of HIV genes (
Figure 1C – black dots), indicates an enrichment of HIV gag, pol and env-expressing cells in latent cells. Merged and individual dominant clusters were analyzed for estimating the effects of DA, alone and in the context of latency. Supplementary Table 1 shows predicted regulators for behaviors in response to latency and DA exposure, as well as genes in each cluster that have been identified in HIV-association studies.

3.1. Uninfected, Unstimulated Control Signatures
Clusters 4, 9, and 11 characterized untreated, uninfected U937 parent promonocytes, referred to as controls (
Figure 1 and
Figure 2). These clusters had a limited number of significant genes analyzed in a rank diagram that indicated overlap in expressed genes, such as the nucleoside diphosphate kinase 2 (NME2) and the ribosomal protein 17 (RPS17), expressed in all 3 clusters (
Figure 2E). The metallothioneins 2A and 1G (MT2A and MT1G), and the amyloid beta precursor protein (APP) were expressed by both dominant clusters 4 and 9, which had a higher number of significant exclusive signatures compared to cluster 11 (
Figure 2D), assigned to pathways important for physiology such as mineral absorption, purine and pyrimidine metabolism, metabolic pathways, neurodegeneration and ribosomes, as seen in circo diagrams (
Figure 2A, 2B and 2C).
3.2. Signatures of DA Stimulation on Uninfected Control Cells
Clusters 1 and 8 were enriched in U937 cells response to DA alone, but not in U1 cells in conditions of HIV latency (
Figure 1B). These clusters overlapped in dimensionality reduction plots, but their profile characterization indicated different phenotypes. Cluster 1 was dominant when compared to cluster 8 due to its higher number of cells, but was phenotypically less active, with a smaller number of genes showing significant expression (
Figure 3A and 3B) compared to cluster 8 (
Figure 3C and 3D). Cluster 8 had 97 significantly upregulated genes and 12 downregulated, while cluster 1 had four out of six significant signatures overlapping with cluster 8 (
Figure 3F). Thus, despite less dominant in numbers, cluster 8 signatures were used to describe effects of DA alone in U937 cells. The four overlapping signature between these clusters were 2 genes downregulated by DA, the Nucleoside Diphosphate Kinase 2 (NME2) and the macrophage migration inhibitory factor (MIF), and 2 upregulated genes, the enzyme cAMP-specific 3',5'-cyclic phosphodiesterase 4D (PDE4B) and the Latent-transforming growth factor beta-binding protein 1 (LTBP1). Two genes were exclusively upregulated in cluster 1, which were SLC1A3 (solute carrier family 1 member 3, glutamate transporter) and PLAC8 (placenta-specific 8). Significant signatures in cluster 8 established networks based on activation, with predicted regulation by TNF and IL4, suggesting activation of inflammatory signals. The signatures in clusters 1 and 8 were linked to pathways in morphine addiction, and metabolism (
Figure 3A and 3C). Common biological processes included response to stimulus, organismal regulation, chemotaxis, hormonal synthesis and response pathways (
Figure 3B and 3D). Cluster 8 was linked to HIV-1 infection, endocytosis and apoptosis pathways, and biological processes linked to signal transduction and immune system-related responses (
Figure 3D). Residual expression of clusters 1 and 8 in U1 cells (
Figure 3G) was correlated with HIV transcriptional activity, as indicated by the dark dots in cluster representations (
Figure 3A and 3D).
3.3. Signatures Restricted to Unstimulated HIV Latently Infected Cells
Clusters 2 and 10 were unique to U1 latent cells and were eliminated by DA stimulation (
Figure 4K), suggesting that these cluster may include critical latency signatures shifted by DA in correlation with HIV transcriptional reactivation (
Figure 4). Cluster 2 was dominant, both in number of events and due to a higher number of significantly expressed genes with significant overrepresentation in pathways with strong gene interactions, indicating orchestrated cellular phenotypes. For instance, cluster 2 phenotypes included genes involved in the integrity of ribosomes (GO:0003735, p=2.8E-23), with 41 out 176 significant signatures that were upregulated ribosome-related structural components (
Figure 4A), predicted to be regulated by the limiting rate eukaryotic initiation factor 5 (EIF5) (
Figure 4B and 4H). Other significant pathways forming strong gene networks were annotated to neurodegeneration (
Figure 4C), metabolic pathways (
Figure 4D), thermogenesis (
Figure 4E), neutrophil extracellular trap formation (
Figure 4F) and oxidative phosphorylation (
Figure 4G), with overlapping and unique signatures and kinases. The neutrophil extracellular trap formation pathway (
Figure 4F) contained integrin ITGAL, voltage-gated channels in mitochondria, and underlying inflammation. These networks indicate the orchestrated behavior of genes interacting via activation, physical binding, reactive and catalytic binding (A, B, R and C connectors, respectively). Molecular processes assigned to cluster 2 indicated that several histone linkers acting on chromatin organization were significantly upregulated. Cluster 10, which also characterized U1 latent cells, had 20 significant genes, 18 downregulated and 2 upregulated, all of them overlapping with cluster 2 (
Figure 4J), but not establishing gene network interactions. Upregulated genes in cluster 10 were the TLR adaptor interacting with endo-lysosomal SLC15A4 (TASL) and proteinase 3 (PRTN3). Moreover, the genes in cluster 10 did not establish any pathway-based interactions.
In addition to EIF5 (
Figure 4B), the three-prime repair exonuclease (MTREX) encoding a protein that removes damaged nucleotides from DNA molecules was also a predicted regulator of latency phenotypes identified in cluster 2 (Supplementary Table 1), while predicted regulators of signatures in cluster 10 were Colony Stimulating Factor 3 (CSF3) and STAT4 (Supplementary Table 1). Both clusters showed representation of signatures that have been observed in gene set representation analyses (GSRA) linked to HIV-1 Gag Pol Nef in immunized humans, and detectable in PBMCs[
43] (M40868, Supplementary Table 1).
3.4. Signatures of DA stimulation on HIV latent cells
The stimulation of U1 latent cells with DA caused enrichment of clusters 3 and 6 (
Figure 5), strongly associated with activation of myeloid functions, including phagocytosis and antigen presentation, as well as pathways linked to viral pathogenesis. Cluster 3 was dominant, due to the higher number of events and due to the higher number of significant signatures annotated to pathways and molecular processes (Figures 5A and 5B, respectively), and establishing pathway-based gene network interactions indicative of orchestrated responses (Figures 5C and 5D). These gene network interactions were annotated to chromosome condensation and remodeling (
Figure 5C) and stress, defense and immune responses (
Figure 5D). Overall, the dominant cluster 3 representing the phenotype of DA-induced latency reversal, was characterized by the downregulation of histone linkers and elements of chromatin condensation and remodeling, as well as components of ribosome structure, assigned to DNA-sensing mechanisms, splicing and pathogenesis. Molecular mechanisms of stress and immune defense included adhesion molecules such as ITGAL and S100A4, and activation markers such as the Allograph Inflammatory Factor 1 (AIF1) also known as Iba1, the Tissue inhibitor of Metalloproteinases 1 (TIMP1) and the apoptosis-associated speck-like protein containing a CARD (PYCARD), indicating a contribution of the inflammasome.
Only 5 signatures were present in cluster 6, all of them overlapped with cluster 3 (
Figure 5G). These overlapping signatures were the upregulation of 2 genes, the long non-coding RNA leucin-rich repeat containing 75A (LRRC75A) and resistin (RETN ) linked to positive regulation of viral replication (
Figure 5F), as well as the downregulation of 3 genes, the FK506 binding protein 5 (FKBP5), the serine-arginine protein kinase 1 (SRPK1) and the brain-expressed X-linked 1 (BEX1), involved in the spliceosome (
Figure 5F).
3.5. Signatures of Latency Reversal Enhanced by DA Stimulation and Their Regulatory Relationships Defined by Gene Network Analysis
The transcription of HIV-1 genes in U1 cells was significantly enhanced by DA exposure (Supplementary Table 2). Computationally, U1 cells with positive transcription of HIV genes Gag-Pol-Env were sorted and their single cell transcriptional profiles were compared to U1 cells negative to HIV genes regardless of DA (
Figure 6).
Signatures of latency reversal included 170 significantly upregulated and 15 downregulated genes. Supplementary Table 2 shows the fifty most upregulated ones, and Supplementary Table 3 shows all the 15 that downregulated ones in cellular events with positive HIV Gag-Pol-Env gene transcription upon comparison to latent cells negative to the expression of the same genes. HIV genes Gag, Pol, Vif, Tat, Env, Vpu, Vpr, Nef and Rev were the most upregulated genes by DA in positive cells (Supplementary Table 2). The examination of the gene list in the Database of Annotation, Visualization and Integrated Discovery (DAVID) resulted in assignments to inflammatory response and immunity (p=1.1E-12), chemotaxis and cell adhesion (p=1E-6), and regulation of transcription (p=1.9E-2) as significant biological processes in phenotypes of active HIV transcription. These phenotypes were also identified in iPathwayGuides, as shown in
Figure 7, with predicted regulators such as TP53, TNF, IL10 and IL4/IL13 involved in both upregulation and downregulation of genes (Figures 7C and 7D). Additional regulators that were identified included lamin B1 (LMNB1) and Sirtuin-1 (SIRT1) (Supplementary Table 1). Supporting the notion that DA is a trigger of the changes linked to latency reversal, neuroactive ligand-receptor interaction was a major pathway identified in these cells, followed by immune and viral pathways (
Figure 7A). Biological processes were centered in immune activation, and protein-containing complex organization, where chromatin structure and spliceosome components are highly represented (
Figure 7B).
For a prediction of HIV latency mechanisms disrupted by DA in the U1 latent model, we compared signatures in dominant clusters characterizing unstimulated latent cells (cluster 2) vs. DA-stimulated latent cells with detectable enhanced viral transcription (cluster 3). The contrasting patterns in gene networks containing ribosomal components of the spliceosome (
Figure 8A and 8B, p=3.9E58), and proteasome (
Figure 8C and 8D, p=4.4E8) visualized by colors, suggested that DA may relax critical mechanisms that explain a decrease in protein-RNA complex organization in protein-containing complexes affecting regulation of inflammation (
Figure 7B) and simultaneously allowing the transcription of viral genes.
The introduction of neighborhood filters using JActiveModules applied to all the identified signatures of latency reversal (170 upregulated genes and 15 downregulated genes, see Supplementary materials) indicated subnetworks with 4 representative patterns (
Figure 9) to indicate contrasting and orchestrated expression that distinguished cluster 2 (dominant latent phenotype) from cluster 3 (dominant latency reversal phenotype), where connectors indicated pathway (blue lines) and physical interactions (red lines).
In Pattern 1, connected genes were significantly downregulated in latent cells (
Figure 9A), but showed contrasting expression levels or no significant changes in clusters cells with positive viral transcription (
Figure 9B). Pattern 1 contained genes involved in biological processes of host-virus interactions (Benjamini p=1.8E-7), apoptosis (p=0.001), and phagocytosis (p=0.009), inclusion of genes with kinase properties (p=8E-4) and common promoter binding motifs to TET2 and SPI1, strongly annotated to pathways of trans-endothelial transmigration (p=5.2E-4).
Pattern 2, on the other hand, was characterized by genes that were significantly overrepresented in latent cells but decreased or no longer significant in cells positive to HIV gene transcription (Figures 9C and 9D). Pattern 2 was linked to biological processes of translation regulation (p=0.05), with molecular functions of ribonucleoprotein (p=9.4E-19) and annotated to the ribosome (p=5E-19) and viral infection pathways (p=1.9E-17).
Pattern 3 was a network of genes downregulated in cells with positive HIV transcription following DA, some of which showed diverse or opposite expression levels in latent cells (Figures 9E and 9F). This subnetwork pattern contained genes also involved in host-virus interactions (p=0.0007), translation regulation (p=0.008) and RNA processing (p=0.02). These genes were linked to molecular functions of helicase (p=0.005), RNA-binding (p=0.005), chaperone (p=0.007) and DNA binding (p=0.05), but links with pathways were not significant. Importantly, patterns 3 and 2 showed substantial overlap in gene representation, which included histone linkers H1-3, H1-4, and H1-5, significantly downregulated in DA-treated cells with positive viral transcription (see also Supplementary Table 2 and Supplementary Figure 2 for protein validation), and with connectors to H1-2 and other histone modifiers in DA-induced pattern 3 (
Figure 9E and 9F). This analysis further suggested the involvement of nuclear components in latency and its reversal caused by DA.
Pattern 4 was characterized by genes mostly significantly upregulated in cells with positive HIV gene transcription, but with inverse or attenuated levels in cells that remain latent (Figures 9G and 9H). The genes in this network are annotated to primary mitochondrial disease (p=0.0004). Biological processes were associated with innate immunity and inflammation (p=0.003), respiratory chain (p=0.008), and molecular functions of oxidoreductase (p=0.04).
3.6. Summary of Results
Overall, the computational strategies and gene network analysis tools applied here to a latent innate immune cell model allowed the identification of interactions linked to latency or its reversal in the context of DA exposure, with annotations that can explain the mechanisms by which DA may affect HIV latent cells in substance users, as summarized in a schematic cartoon in
Figure 10.
4. Discussion
The goals of this study were to enhance our understanding of HIV latency relative to active transcription in myeloid cells and to explore the effects of external factors, particularly dopamine (DA), on transcriptional regulation. Reports of increased viral load in the brains of stimulant substance users suggest that HIV latency mechanisms can be disrupted by complex factors in these individuals [
44]. Given the well-known relationship between psychostimulants and DA, our experimental design examined the transcriptional profiles of control and latently infected promonocytes to uncover DA-induced gene expression changes that could explain latency reversal in substance users. We previously demonstrated that DA promotes HIV entry into innate immune cells [
1] and enhances viral transcription in a latency model [
2]. However, DA activates HIV gene transcription in only a subset of cells. Single-cell RNA sequencing enabled us to isolate these responsive subpopulations and compare them to uninfected, latently infected, and DA-stimulated cells within a promonocyte cell line system. This approach lays the groundwork for future in vivo studies.
A series of analytical approaches proved useful for identifying cellular clusters expressing specific behaviors as a result of DA stimulation and in the context of HIV. A strength of the introduction of systems strategies is the detection of biological processes rather than single gene signatures, increasing the power of discovery of redundancies when response patterns are characterized by more than one cluster. For instance, clusters 4, 9, and 11 characterized unstimulated U937; clusters 1 and 13 characterized U937 cells stimulated with DA; clusters 2 and 10 characterized unstimulated U1 latent cells; clusters 3, 6, and 12 were found in U1 latent cells stimulated with DA. The separation in clusters results from differences in signatures ranked by dominance. Although common signatures linked to the overall cluster behaviors could be identified in some cases, systems interactions and annotated pathway and biological process differences were more useful to understand the role of DA as a factor reversing latency.
The significance of our approach is enhanced by the lack of overall understanding of cellular factors and cellular pathways mediating HIV latency and reactivation. In the U1 cell model used in this study, an attenuated Tat confers the characteristic of latency [
29,
45]. Yet, PKC activators and latency reversal agents strongly increase p24 in these cells’ supernatant because of intact and viable HIV replication in their totality [
2,
46], which ensures that the model is valid for developing a study pipeline. Indeed, we were able to detect patterns that characterize latency and HIV transcription in the context of DA exposure, with signatures that provided mechanistic insights.
The uninfected parent promonocytic cells displayed a limited number of distinctive transcriptional signatures, notably increased expression of nucleoside diphosphate kinase B (NME2), which facilitates the synthesis of nucleoside triphosphates other than ATP, and Ribosomal Protein S17 (RPS17), as well as upregulation of amyloid beta precursor protein (APP) and the metallothioneins MT1G and MT2A. These cysteine-rich metallothioneins play key roles in DNA damage repair through metal detoxification, modulation of mineral absorption, and enhancement of cellular resistance to apoptosis. Consistent with our observations, previous studies have highlighted the essential function of these molecules in preserving normal monocyte activity [
47]. NME2 is a gene encoding a protein complex that helps with trafficking between the endoplasmic reticulum and Golgi and that may have epigenetic functions depending on its tertiary structure[
48]. DA stimulation on uninfected conditions on the other hand, caused downregulation of NME2 and MIF and the upregulation of PDE4B and LTBP1. The contrasting effect on NME2 expression in control cells following DA indicated a potential signature of DA stimulation in promonocytes.
The dominant phenotype in unstimulated U1 cells characterizing latency involved the high expression of genes encoding ribosome structural elements, as well as replication-dependent histones such as H2AC17 and histone linkers. These findings indicated the contribution of chromatin structure in transcriptional silencing, while indicating that protein-containing complex organization, with a predicted regulatory role of the eukaryotic translation initiation element EIF5, are linked to HIV latency. Interestingly, the signatures of latency were assigned to neurodegeneration, metabolic pathways, oxidative phosphorylation networks, and the presence of molecules controlling leukocyte transmigration such as ITGAL, which indicate despite the latent phenotype, these cells may have the capacity to influence the environment to enhance inflammation. This could be a result of residual viral replication[
49] and could explain the persistence of low levels of neuroinflammation in ART-suppressed PWH.
The effects of DA on HIV latent cells were detectable in the cluster 3, where myeloid and inflammatory pathways were significantly enhanced, including PYCARD, an adaptor component of the inflammasome previously reported to be activated by DA in monocytes[
50]. Importantly, genes encoding histone linkers and chromatin structural elements were significantly downregulated by DA in latent cells. The computational sorting of cells positively transcribing HIV genes indicated that a decrease in histone linkers is a strong signature of DA-induced latency reversal, suggesting effects on chromatin organization. Moreover, genes linked to stress response, and immune defense, including anti-viral response and inflammation, were enhanced in these cells. The predicted regulators of latency reversal phenotypes by DA included TP53, TNF, as well as IL10 and the IL13/IL4 axis. IL13 and IL4 can control myeloid phenotypes towards suppressor M2 responses [
51,
52,
53]. Overall, DA may shift the response of latent myeloid cells towards activation, with the contribution of changes of chromatin structure, due to the decrease in histone linkers H1-5, H1-4, opening chromatin (see Supplementary materials) and unrepressing HIV transcription along with pro-inflammatory genes. While others have suggested that DA is a suppressor of inflammation[
52,
53], our previous studies in other innate immune cell models in the context of HIV infection, as well as in U1 cells containing latent virus, have indicated that DA can activate a significant subset of cells in a dopamine receptor signaling-dependent manner [
1].
An additional strategy involved gene network analysis, for identifying genes with relationships based on pathway and physical interactions, and showing orchestrated behaviors, and expression patterns that could distinguish between latent and reversed HIV transcription in the context of DA exposure. Overexpression of ribosomal and chromatin structure components remained as signatures of interest in latent cells, while latency reversal was linked to inflammation and mitochondrial respiratory chain. The gene networks with connectors based exclusively on pathways were annotated to ribosomal control of transcription and the proteasome, overall assigned to protein-containing complex organization. When physical interactions were added, biological processes were expanded to viral-host interactions, mitochondrial health, and viral infection, supporting the notion that cells with transcriptionally active HIV expression have an active phenotype[
54]. Genes annotated to dopaminergic signaling were also identified as perturbed by DA in correlation with increased HIV transcription, further indicating specific DA signaling involved in the activation of these cells.
Additional considerations may be introduced regarding the effects of DA on epigenetic silencers such as SIRT1, which was predicted as a regulator of DA-induced phenotypes, and play a role in the stabilization of histone linkers. We previously described a loss of Sirt1 and its epigenetic silencing properties as a cause of neuroinflammation in microglia using the non-human primate model of neuroHIV[
55]. The results presented here suggest that chromatin structural factors may underlie reactivation and may be considered for intervention approaches.
While the U1 promonocytic cell line is a well-established model for studying HIV latency, it is important to acknowledge the intrinsic limitations associated with the use of immortalized cell lines. The U1 cells, derived from the U937 monocyte lineage, may not fully represent the phenotypic diversity found in primary myeloid cells within the central nervous system of HIV-infected individuals. Furthermore, the cell culture environment lacks the complexity of tissue architecture and systemic influences present in vivo. However, despite the identification of two mutations in Tat causing a partial deficiency that is regarded as a driver of latency[
56], baseline p24 is detectable (Supplementary Figure 1) and U1 cells can be induced to produce virus, which in this study was confirmed in the supernatants of the analyzed cultures by p24 ELISA, and further validating the model for studying HIV latency and reactivation. Apoptotic pathways mediated by p53 expression in U1 cells and ectopic Tat can also reactivate HIV-1 transcription and replication from its latent state[
57]. However, this model has been also used to examine mechanisms of Tat-independent transcription [
57]. Similarly to cells in vivo, cytokines such as interferon-gamma (IFN-γ)[
58] also can activate HIV expression in U1 cells, including redirecting the production of virions to intracytoplasmic vacuoles [
58]. The HIV produced by U1 cells cannot infected new cells despite producing p24[
59]. Yet, whether the resulting virus is infectious or not does not decrease the value of the model for predicting mechanisms of reactivation, and the phenotypes that allow the identification of latent states. Other models should be used to investigate the effects of DA on infection, as previously published by us[
1,
2]. U1 cells as a HIV latent model are useful to our approach aiming at making predictions of mechanisms and especially unique signatures to help identify latent cells in systems where the virus is not tagged by a reporter or within the heterogeneity of target cells in vivo, and search best strategies to reverse latency by harnessing the immune-dopaminergic response. The effects of DA enhancing viral transcription are replicated in other cellular models of HIV latency that are under investigation, including in primary human microglia and in C20 microglia cell line (Supplementary Figure 1). Most importantly, the model provided the development of an analytical strategy guideline that is applicable to other models. The mechanistic insights and the method have utility in validation, confirmation and for expanding the knowledge gained here. We have confirmed that DA activates myeloid cells.
The single cell approach allowed us to observe that the effects are heterogeneous, even in cell lines, with a significant subset of responding cells, detectable in clusters with viral transcription. U1 cells may be intrinsically heterogeneous due to the replication competence of the integrated virus causing dynamism in integration sites[
60,
61], which could potentially explain the phenotypic clusters and response pattern diversity by condition. However, heterogeneity in cell lines has been previously detected using single cell strategies [
62]. Moreover, uninfected parent U937 cells also showed heterogeneity, indicating that the proviral integration sites are not solely responsible for clusters defined by differences in transcriptional signatures. It is important to acknowledge that there are no perfect models, but the within-cell line diversity is a strength, for emulating the diversity of HIV innate immune target cells in the brain, including microglia[
63,
64,
65]. Future studies employing primary human myeloid cells or animal models are warranted to validate the findings reported here, especially applicable to microglia, and to ascertain their relevance in a physiological context.
Whether the reactivation of HIV genes is part of a general effect on gene transcription leading to inflammation caused by DA, but that may occur through other means regardless of DA, remains to be confirmed. The decrease in histone linkers in correlation with active transcription of HIV genes indicates a critical role for DA in nucleosome components and organization, underlying changes in the epigenetic landscape in the context of substance use, as described by us and others[
55,
66,
67]. Overall, the strategy allowed the identification of a subset of cells responding to DA, and their signatures linked to HIV gene transcription, which can be applied to more complex model systems and brains from PWH that are substance users.
Figure 1.
– UMAP plots of transcriptomic data in single cell RNAseq datasets. U937 and U1 cells were cultured with vehicle or 1uM DA for 24 hrs. Cells were collected and processed for characterization using single-cell RNAseq. Cluster analysis and visualization were performed. (A) UMAP for dimention reduction of data in individual treatment conditions and (B) merged treatment conditions with color-coded clusters. (C) Individual clusters with dark dots indicating events positive to HIV Gag-Pol-Env genes in U1 cell populations.
Figure 1.
– UMAP plots of transcriptomic data in single cell RNAseq datasets. U937 and U1 cells were cultured with vehicle or 1uM DA for 24 hrs. Cells were collected and processed for characterization using single-cell RNAseq. Cluster analysis and visualization were performed. (A) UMAP for dimention reduction of data in individual treatment conditions and (B) merged treatment conditions with color-coded clusters. (C) Individual clusters with dark dots indicating events positive to HIV Gag-Pol-Env genes in U1 cell populations.
Figure 2.
- Clusters 4, 9 and 11 characterizations, overlapping signatures, and rank diagram. A) Circos plot showing significant genes in cluster 4 as per Bonferroni adjusted p values and pathway assignment visualization. B) Circos plot showing significant genes in cluster 9 as per Bonferroni adjusted p values and pathway assignment visualization. C) Circos plot showing significant genes in cluster 11 as per Bonferroni adjusted p values and pathway assignment visualization. D) Gene intersection analysis between clusters 4, 9 and 11, E) Venn’s Diagram indicating signature overlap.
Figure 2.
- Clusters 4, 9 and 11 characterizations, overlapping signatures, and rank diagram. A) Circos plot showing significant genes in cluster 4 as per Bonferroni adjusted p values and pathway assignment visualization. B) Circos plot showing significant genes in cluster 9 as per Bonferroni adjusted p values and pathway assignment visualization. C) Circos plot showing significant genes in cluster 11 as per Bonferroni adjusted p values and pathway assignment visualization. D) Gene intersection analysis between clusters 4, 9 and 11, E) Venn’s Diagram indicating signature overlap.
Figure 3.
- Clusters 1 and 8 characterization, pathway and molecular process assignments, and linked gene network analysis. A) Window showing Cluster 1 and a Circos plot showing significant genes in cluster 1 as per Bonferroni adjusted p values and pathway assignment visualization. B) Circos plot showing significant molecular processes in cluster 1. C) Circos plot showing significant genes in cluster 8 as per Bonferroni adjusted p values and pathway assignment visualization. D) Window showing Cluster 8 and Circos plot showing significant biological processes. E) Networks established by genes assigned to significant pathways in cluster 8, with connectors indicating interactions based on activation, and predicted regulators. Pink nodes indicate upregulated genes, and blue nodes indicate downregulated genes. F) Cluster 1 and 8 directional signature overlap. G) Cluster 1 and 8 were enriched by DA alone in U937 uninfected parent cells, but occurred is small numbers in U1 latent cells.
Figure 3.
- Clusters 1 and 8 characterization, pathway and molecular process assignments, and linked gene network analysis. A) Window showing Cluster 1 and a Circos plot showing significant genes in cluster 1 as per Bonferroni adjusted p values and pathway assignment visualization. B) Circos plot showing significant molecular processes in cluster 1. C) Circos plot showing significant genes in cluster 8 as per Bonferroni adjusted p values and pathway assignment visualization. D) Window showing Cluster 8 and Circos plot showing significant biological processes. E) Networks established by genes assigned to significant pathways in cluster 8, with connectors indicating interactions based on activation, and predicted regulators. Pink nodes indicate upregulated genes, and blue nodes indicate downregulated genes. F) Cluster 1 and 8 directional signature overlap. G) Cluster 1 and 8 were enriched by DA alone in U937 uninfected parent cells, but occurred is small numbers in U1 latent cells.
Figure 4.
– Cluster 2 characterization, assignments, and gene network analysis. A) Circos plot showing significant genes in cluster 2 as per Bonferroni adjusted p values and pathway assignment visualization. B, C, D, E, F, G) Networks established by genes assigned to each significant pathway in cluster 2, with connectors indicating interactions based on activation, binding, catalysis, inhibition, and reaction. Pink nodes indicate genes with increased fold change (FC), and blue nodes indicate decreased FC in relation to average expression. B) Ribosome pathway gene interactions, C) Neurodegeneration, D) Metabolic Pathways, E) Thermogenesis, F) Neutrophil extracellular trap formation pathway, G) Oxidative Phosphorylation. H) Circos plot showing significant molecular processes in cluster 2. I) Gene network analysis for non-pathway overlapping genes assigned to structural molecular activity, and nucleic acid binding. J) Clusters 2 and 10 signature overlap, and representation of dark dots indicating HIV Gag-Pol-Env-positive events occurring in U1 cells, and their directional phenotypic intersect, indicating dominant cluster 2 patterns. K) Map of merged clusters 2 and 10 unique to U1 latent cells and their pattern suppression in response to DA.
Figure 4.
– Cluster 2 characterization, assignments, and gene network analysis. A) Circos plot showing significant genes in cluster 2 as per Bonferroni adjusted p values and pathway assignment visualization. B, C, D, E, F, G) Networks established by genes assigned to each significant pathway in cluster 2, with connectors indicating interactions based on activation, binding, catalysis, inhibition, and reaction. Pink nodes indicate genes with increased fold change (FC), and blue nodes indicate decreased FC in relation to average expression. B) Ribosome pathway gene interactions, C) Neurodegeneration, D) Metabolic Pathways, E) Thermogenesis, F) Neutrophil extracellular trap formation pathway, G) Oxidative Phosphorylation. H) Circos plot showing significant molecular processes in cluster 2. I) Gene network analysis for non-pathway overlapping genes assigned to structural molecular activity, and nucleic acid binding. J) Clusters 2 and 10 signature overlap, and representation of dark dots indicating HIV Gag-Pol-Env-positive events occurring in U1 cells, and their directional phenotypic intersect, indicating dominant cluster 2 patterns. K) Map of merged clusters 2 and 10 unique to U1 latent cells and their pattern suppression in response to DA.

Figure 5.
- Cluster 3 and 6 pathway and molecular process assignments and linked gene network analysis. A) Circos plot showing significant genes in cluster 3 as per Bonferroni adjusted p values, and pathway assignment visualization. B) Circos plot showing biological processes assignments in cluster 3. C) Gene network containing genes overlapping between chromatin condensation and remodeling biological processes. D) Gene network containing overlapping genes in response to stress, defense and immune response processes. Connectors indicate interactions based on activation, binding, catalysis, inhibition and reaction. Pink nodes indicate genes with increased FC, and blue nodes indicate decreased FC in relation to average expression. E) Circos plot showing Cluster 6 pathway assignments for significant genes. F) Circos plot showing Cluster 6 Biological processes. G) Clusters 3 and 6 directional intersect, indicating dominant cluster 3 pattern. K) Clusters 3 and 6 distributions unique to U1 latent cells stimulated by DA.
Figure 5.
- Cluster 3 and 6 pathway and molecular process assignments and linked gene network analysis. A) Circos plot showing significant genes in cluster 3 as per Bonferroni adjusted p values, and pathway assignment visualization. B) Circos plot showing biological processes assignments in cluster 3. C) Gene network containing genes overlapping between chromatin condensation and remodeling biological processes. D) Gene network containing overlapping genes in response to stress, defense and immune response processes. Connectors indicate interactions based on activation, binding, catalysis, inhibition and reaction. Pink nodes indicate genes with increased FC, and blue nodes indicate decreased FC in relation to average expression. E) Circos plot showing Cluster 6 pathway assignments for significant genes. F) Circos plot showing Cluster 6 Biological processes. G) Clusters 3 and 6 directional intersect, indicating dominant cluster 3 pattern. K) Clusters 3 and 6 distributions unique to U1 latent cells stimulated by DA.

Figure 6.
– Distribution of HIV gene transcripts in U1 cells at baseline and following DA stimulation for 24 hrs. The colored scale bar denotes the relative expression values of HIV-1 genes Gag-Pol (left panels) and Env (right panels), at vehicle baseline conditions (upper panels) and after the treatment with DA (lower panels).
Figure 6.
– Distribution of HIV gene transcripts in U1 cells at baseline and following DA stimulation for 24 hrs. The colored scale bar denotes the relative expression values of HIV-1 genes Gag-Pol (left panels) and Env (right panels), at vehicle baseline conditions (upper panels) and after the treatment with DA (lower panels).
Figure 7.
– Pathway and biological processes annotations, and predicted regulators to signatures that differ between U1 cells with active replication following DA stimulation vs. cells that preserve latency. The comparison between cells transcribing HIV genes and cells that do not resulted in significantly increased and decreased gene signatures. A) Circos plot indicating main annotated pathways. B) Circos plot indicating main biological processes linked to gene signatures. D) Predicted regulators of upregulated genes. E) Predicted regulators of downregulated genes. E-H) Gene network interactions with consensus predicted regulators, where white nodes represent predicted regulators, pink nodes represent genes upregulated and blue downregulated in correlation with active HIV transcription. E) TP53-regulated gene signatures, F) TNF-regulated genes, G) IL10-regulated genes, H) IL4/IL13 regulated genes.
Figure 7.
– Pathway and biological processes annotations, and predicted regulators to signatures that differ between U1 cells with active replication following DA stimulation vs. cells that preserve latency. The comparison between cells transcribing HIV genes and cells that do not resulted in significantly increased and decreased gene signatures. A) Circos plot indicating main annotated pathways. B) Circos plot indicating main biological processes linked to gene signatures. D) Predicted regulators of upregulated genes. E) Predicted regulators of downregulated genes. E-H) Gene network interactions with consensus predicted regulators, where white nodes represent predicted regulators, pink nodes represent genes upregulated and blue downregulated in correlation with active HIV transcription. E) TP53-regulated gene signatures, F) TNF-regulated genes, G) IL10-regulated genes, H) IL4/IL13 regulated genes.
Figure 8.
– Pathway-based connections among gene signatures significantly perturbed by DA in U1 latent cells. Pathway-based gene networks identified in (A, C) clusters 2 and 10 that characterize U1 latent cells, in comparison to (B, D) clusters 3 and 6 that characterized DA-treated cells exhibiting latency reversal. The most significant pathway annotations were (A, B) ribosomal control of transcription (Benjamini, p=3.9E-58) and (C, D) proteasome (Benjamini, p=4.4E-8). Colors indicate upregulation (shades of red) and downregulation (shades of blue). Blue connectors indicate links through pathway.
Figure 8.
– Pathway-based connections among gene signatures significantly perturbed by DA in U1 latent cells. Pathway-based gene networks identified in (A, C) clusters 2 and 10 that characterize U1 latent cells, in comparison to (B, D) clusters 3 and 6 that characterized DA-treated cells exhibiting latency reversal. The most significant pathway annotations were (A, B) ribosomal control of transcription (Benjamini, p=3.9E-58) and (C, D) proteasome (Benjamini, p=4.4E-8). Colors indicate upregulation (shades of red) and downregulation (shades of blue). Blue connectors indicate links through pathway.
Figure 9.
– Subnetworks with orchestrated expression patterns and contrasting differences between cluster sets representing latent cells in resting state and reversed latency following DA stimulation. Subnetworks were identified in the merged pathway and physical interactions by applying filters using JActiveModules in Cytoscape 3.10.2. Subnetworks indicated 4 expression patterns that allowed a parallel visualization of contrasting behaviors between dominant response clusters. Cluster 2 representing the dominant latency phenotype (A, C, E, G) was visualized in parallel to cluster 3 (B, D, F, H) representing the dominant DA-induced latency reversal phenotype. The identified modules reflected (A, B) Pattern 1 containing 72 nodes, (C, D) Pattern 2 containing 21 nodes, (E, F) Pattern 3 containing 48 nodes, and (G, H) Pattern 4 with 59 nodes. Red connector lines indicate physical interactions. Blue connectors indicate pathway interactions. Transcriptional expression levels are indicated by node colors, with shades of blue indicating downregulated and shades of red indicating upregulated genes. Diamond shapes indicate statistical significance (adjusted p<0.05), and circles indicate non-significant genes.
Figure 9.
– Subnetworks with orchestrated expression patterns and contrasting differences between cluster sets representing latent cells in resting state and reversed latency following DA stimulation. Subnetworks were identified in the merged pathway and physical interactions by applying filters using JActiveModules in Cytoscape 3.10.2. Subnetworks indicated 4 expression patterns that allowed a parallel visualization of contrasting behaviors between dominant response clusters. Cluster 2 representing the dominant latency phenotype (A, C, E, G) was visualized in parallel to cluster 3 (B, D, F, H) representing the dominant DA-induced latency reversal phenotype. The identified modules reflected (A, B) Pattern 1 containing 72 nodes, (C, D) Pattern 2 containing 21 nodes, (E, F) Pattern 3 containing 48 nodes, and (G, H) Pattern 4 with 59 nodes. Red connector lines indicate physical interactions. Blue connectors indicate pathway interactions. Transcriptional expression levels are indicated by node colors, with shades of blue indicating downregulated and shades of red indicating upregulated genes. Diamond shapes indicate statistical significance (adjusted p<0.05), and circles indicate non-significant genes.

Figure 10.
– Summary of findings - Innate immune cells are maintained in homeostatic resting state with the expression of strong signatures such as NME2 and MIF, and keeping physiological processes via a healthy purine metabolism and mineral absorption (upper left). In cells that bear HIV provirus and are virologically suppressed, silence is maintained by a high expression of histone linkers and other components of chromatin organization, at the cost of a higher metabolic activation that may lead to some neurodegeration (lower left). In the context of high levels of DA, innate immune cells do increase pathways associated with addiction, mainly dopaminergic signaling due to the presence of DA receptors in the innate immune cell surface, and correlated with a decrease of the latency signature genes NME2 and MIF (upper right). On the other hand, cells with a provirus that are exposed to DA are strongly affected with a significant decrease in histone linkers and protein-containing complex organization, unrepressing inflammation and unleashing HIV transcription. Predicted regulators of these phenotypes include EIF5, TP53, immune cytokines such as TNF, IL10, IL4 and IL13, as well as epigenetic silencers such as SIRT1.
Figure 10.
– Summary of findings - Innate immune cells are maintained in homeostatic resting state with the expression of strong signatures such as NME2 and MIF, and keeping physiological processes via a healthy purine metabolism and mineral absorption (upper left). In cells that bear HIV provirus and are virologically suppressed, silence is maintained by a high expression of histone linkers and other components of chromatin organization, at the cost of a higher metabolic activation that may lead to some neurodegeration (lower left). In the context of high levels of DA, innate immune cells do increase pathways associated with addiction, mainly dopaminergic signaling due to the presence of DA receptors in the innate immune cell surface, and correlated with a decrease of the latency signature genes NME2 and MIF (upper right). On the other hand, cells with a provirus that are exposed to DA are strongly affected with a significant decrease in histone linkers and protein-containing complex organization, unrepressing inflammation and unleashing HIV transcription. Predicted regulators of these phenotypes include EIF5, TP53, immune cytokines such as TNF, IL10, IL4 and IL13, as well as epigenetic silencers such as SIRT1.
