Histone deacetylase inhibitors induce complex host responses that contribute to differential potencies of these compounds in HIV reactivation

Histone deacetylase (HDAC) inhibitors (HDACis) have been widely tested in clinical trials for their ability to reverse HIV latency but have yielded only limited success. One HDACi, suberoylanilide hydroxamic acid (SAHA), exhibits off-target effects on host gene expression predicted to interfere with induction of HIV transcription. Romidepsin (RMD) has higher potency and specificity for class I HDACs implicated in maintaining HIV provirus in the latent state. More robust HIV reactivation has indeed been achieved with RMD use ex vivo than with SAHA; however, reduction of viral reservoir size has not been observed in clinical trials. Therefore, using RNA-Seq, we sought to compare the effects of SAHA and RMD on gene expression in primary CD4+ T cells. Among the genes whose expression was modulated by both HDACi agents, we identified genes previously implicated in HIV latency. Two genes, SMARCB1 and PARP1, whose modulation by SAHA and RMD is predicted to inhibit HIV reactivation, were evaluated in the major maturation subsets of CD4+ T cells and were consistently either up- or down-regulated by both HDACi compounds. Our results indicate that despite having different potencies and HDAC specificities, SAHA and RMD modulate an overlapping set of genes, implicated in HIV latency regulation. Some of these genes merit exploration as additional targets to improve the therapeutic outcomes of “shock and kill” strategies. The overall complexity of HDACi-induced responses among host genes with predicted stimulatory or inhibitory effects on HIV expression likely contributes to differential HDACi potencies and dictates the outcome of HIV reactivation.

ities, SAHA and RMD modulate an overlapping set of genes, implicated in HIV latency regulation. Some of these genes merit exploration as additional targets to improve the therapeutic outcomes of "shock and kill" strategies. The overall complexity of HDACi-induced responses among host genes with predicted stimulatory or inhibitory effects on HIV expression likely contributes to differential HDACi potencies and dictates the outcome of HIV reactivation.
Clinical trials using histone deacetylase (HDAC) 4 inhibitors (HDACis) to reverse HIV latency have shown limited success in demonstrating induction of HIV gene expression in CD4 ϩ T cells of individuals on suppressive anti-retroviral therapy (ART). Following treatment of cells from HIV-infected individuals ex vivo with suberoylanilide hydroxamic acid (SAHA), induction of HIV transcription was demonstrated initially (1)(2)(3)(4), and subsequently induction of HIV protein production was demonstrated both ex vivo (5) and following in vivo administration (6). Although cells with HIV protein expression induced by SAHA were killed by cytotoxic T lymphocytes ex vivo (5,6), reduction in the size of the latent reservoir following administration of SAHA in vivo was not observed (2,3). Another HDACi, romidepsin (RMD) appeared to be more potent than SAHA for reactivating HIV ex vivo (7); however, administration of RMD alone in vivo did not result in a reduction of the latent reservoir (4). A more recent clinical trial using a combined treatment strategy of RMD with Vacc-4x, recombinant human granulocyte macrophage colony-stimulating factor vaccination, reported a reduction of the reservoir 6 weeks after RMD administration (8). This reduction in reservoir size was, however, not statistically significant (8).
By interfering with the active site of HDAC enzymes, HDA-Cis cause hyperacetylation of histones and chromatin relaxation. HDAC classification into four distinct classes is based on homology to yeast acetylases. HDACs from classes I, II, and IV are Zn 2ϩ -dependent enzymes, whereas class III HDACs, called sirtuins, are NAD ϩ -dependent (reviewed in Ref. 9). The HDA-Cis, used in "shock and kill" strategies to reactivate latent HIV, inhibit HDACs of classes I, II, and IV by chelating Zn 2ϩ ions or becoming covalently linked at the active site (reviewed in Ref. 9), with different potencies and specificities (7,10). A number of HDACis have a broad range of action (called pan-HDACis); however, the IC 50 values of different HDACs vary. For example, SAHA is active against representative HDACs from classes I, II, and IV; however, it is most potent against HDAC6, a representative of class II (7). In contrast, RMD has high potency against HDACs of class I with very low activity against HDAC6 (7,11). Because HDACs of class I have been implicated in maintaining HIV provirus in a latent state (12)(13)(14), higher potency of RMD for these particular HDACs may explain why it was more active than SAHA for reactivation of HIV.
Despite chromatin relaxation caused by HDACi, gene expression-profiling studies of cells treated with SAHA identified comparable numbers of up-and down-regulated genes (15)(16)(17). These findings are consistent with the idea that alternative mechanisms of gene regulation are induced by SAHA, such as secondary down-modulation of gene expression by activated transcription factors or modulation mediated via nonhistone effects (reviewed in Ref. 18). A recent study by Zaikos et al. (19) demonstrated that SAHA induced inhibitory effects on HIV reactivation, specifically via its activity against HDAC6 and its nonhistone targets. RMD exhibited a similar inhibitory behavior on HIV expression (19), although it has only marginal activity against HDAC6 (7,11). It is therefore plausible that the inhibitory effect of RMD on HIV is not associated with direct inhibition of HDAC6 activity. It is likely that secondary effects on host gene expression following the initial gene activation play an important role in the net achievable level of induction of HIV gene expression following treatment with HDACi. We have previously shown that some of these secondary effects in the case of SAHA represented down-regulation of HIV transcriptional activators and up-regulation of repressors, both at the RNA and protein levels (20). The present study was undertaken to assess the effect of RMD on the transcriptome of primary uninfected CD4 ϩ T cells and an in vitro model of HIV latency and to compare the effects of RMD and SAHA on host genes to enhance our understanding of the secondary mechanism of action of these compounds that are relevant to reactivation of latent HIV. A subset of host genes modulated by these two HDACis were assessed as potential targets for improving shock and kill treatment outcomes.

Choice of HDACi doses and verification of activity
Treatment with 1 M SAHA has been shown previously to result in induction of expression of HIV RNA in the "Spina model" of HIV latency (21) that was used in the present study. Here, we have tested the activity of RMD for HIV reactivation in the same model system. To determine the optimal concentration of RMD, cells from the latency model were treated with different concentrations of RMD, and HIV expression (gag RNA) was measured by droplet digital PCR (ddPCR). We found that RMD had similar activity around a range of doses, with the average peak of activation observed at 15 nM (Fig. 1A). Immunoblotting using an antibody against acetylated histone H3 demonstrated that at selected doses, treatment with SAHA or RMD resulted in an increased acetylation of histone H3 (Fig.  1B), confirming that both compounds were active.

Viability, T cell activation, and induction of HIV transcription following treatment with HDACi
Paired samples of the model of HIV latency and mock-infected cells, treated with 1 M SAHA, 15 nM RMD, or their solvent control DMSO for 24 h, were collected for RNA-Seq. With selected treatment conditions, treatment with SAHA or RMD did not reduce cell viability, compared with DMSO treatment ( Fig. 2A) as assessed by trypan blue staining. Cell activation was assessed by measuring surface proteins expressed during early (CD69), intermediate (CD25), and late (HLA-DR) stages of T cell activation; cell exhaustion was assessed using surface expression of CD279 (also known as PD-1). In general, the percentage of cells positive for any of these markers, following treatment, was not greater than that for the untreated cells.
In the in vitro model of latency, there were significantly more CD69-positive cells in RMD-treated samples, compared with DMSO-treated (average -fold increase 1.6, p ϭ 0.0006). However, this difference appeared to be caused by the exposure to DMSO, which resulted in reduction in frequency of CD69-positive cells, compared with untreated sample (average -fold decrease 1.6, p Ͻ 0.0001). On the other hand, exposure of cells to virus during establishment of the in vitro model significantly reduced viability (p Ͻ 0.0001) and induced modest T cell acti- A, dose-response curve of RMD treatment of the in vitro model of HIV latency to determine the optimal dose for HIV reactivation. The experiment was performed using cells from five independent donors. Data are presented as mean -fold changes of HIV gag expression (RMD/ DMSO). Error bars, S.D. B, immunoblot using an antibody against acetylated histone H3 (Ac-H3) of cell extracts treated with DMSO, 1 M SAHA, or 15 nM RMD for 24 h. Lamin B1 was measured as a sample-loading control. The experiment was performed using uninfected cells from three independent donors.
HIV transcription was variably induced in the four sequenced samples, as assessed by ddPCR. Induction of gag RNA ranged from no response in one donor to a 1.85-fold increase with SAHA and a 2.75-fold increase with RMD (Fig. 2B). Similar trends were seen with total polyadenylated (poly(A)) and singly spliced (env) transcripts. Multiply spliced (MS) transcripts were the least responsive in this experiment (Fig. 2B). The degree of HIV reactivation did not appear to correlate with the percentage of latently infected cells. Of the two samples with the most integrants (in 20% of cells), one was the best responder to the HDACi, whereas the other one was the nonresponder (Fig. 2C).

Exposure of cells to virus does not significantly alter the response of cells to HDACi
Because the frequency of latently infected cells is low in our in vitro model (8 -20%), any difference in their response to HDACi would unlikely be detectable, unless it was very robust (i.e. Ͼ5-10-fold). However, it is not known whether exposure of cells to virus during the culture period affects the way they respond to HDACi treatment. Because exposure to virus slightly affected the activation status of the cells, we sought to test whether cells in the model of HIV latency differ in response to HDACi, compared with mock-infected cells. A complex contrast approach was used in EdgeR to identify genes, differentially modulated by the two HDACis (see supporting Methods for details). This formal analysis demonstrated that only 6 and 11 genes were differentially modulated by SAHA and RMD, respectively, when responses in mock-infected cells were compared with cells from the HIV latency model (Table S1). These results are consistent with the idea that the effect of exposure to virus on alterations of gene expression induced by HDACi is minimal and that gene modulations were similar in mock-infected cells and the model of HIV latency. One of the differentially modulated genes was cyclin-dependent kinase inhibitor 1A (CDKN1A/p21), a well-known host factor induced by many HDACi, including SAHA and RMD. Its up-regulation by HDACi mediates growth inhibition and apoptosis in cancer cells (22)(23)(24). In our study, CDKN1A was highly (-fold increase

Secondary mechanisms of action of HDACi and HIV latency
of 3.6 for SAHA and 3.9 for RMD) and significantly up-regulated in mock-infected cells (false discovery rate corrected p value (FDR) Ͻ 0.05) but only minimally up-regulated (-fold increase of 1.28 and FDR Ͼ 0.05) in cells from the HIV latency model (Table S1).

Host genes and cellular pathways affected by SAHA and RMD
We next sought to identify human genes modulated by SAHA and RMD in primary CD4 ϩ T cells to better understand what effects on the host may be stimulatory or inhibitory for HIV reactivation. We previously performed a gene expression profiling study, using microarrays coupled with LC-MS to identify genes modulated by SAHA (20). Here, we took advantage of the sensitivity of RNA-Seq technology to detect a greater number of genes not identified previously (Table S1 shows all genes that were detected in this study, with -fold changes and p values). To provide a comparison, this study found 11,596 genes significantly modulated by SAHA (FDR Ͻ 0.05) in mock-infected cells of 16,059 total genes detected, whereas the microarray study identified 2,982 genes (20). Because many more genes were identified, for the purpose of the subsequent analyses, we have filtered genes based on both significance (FDR Ͻ 0.05) and absolute -fold change (Ͼ1.5). Under these filtering criteria, in mock-infected cells, 3,440 genes (34.5%) were up-regulated and 2,550 genes (25.6%) were down-regulated in common by both HDACis. In the model of latency, 3,089 genes (31.1%) were up-regulated and 2,506 genes (25.3%) were down-regulated in common by both HDACis. Of the genes modulated in common, the magnitude of response was generally greater in cells treated with RMD than with SAHA (72.3-82.5% of genes had greater responses to RMD, up or down, both in mock-infected cells and in the model of HIV latency). RMD also induced expression changes in a larger number of additional unique genes than did SAHA. For example, with cells from the latency model, RMD up-regulated (FDR Ͻ 0.05, -fold change Ͼ 1.5) 1,975 genes (19.9%) that did not respond to SAHA, whereas only 405 genes (4.1%) were up-regulated by SAHA and not RMD (Fig. 3). These relationships between treatment responses were similar for the down-regulated genes in cells of the latency model and for all genes modulated in mock-infected cells (Fig.  3). A small fraction of genes appeared to be modulated in oppo-site directions by SAHA and RMD; 54 genes (0.5%) were upregulated by RMD and down-regulated by SAHA, and 5 genes (0.1%) were up-regulated by SAHA and down-regulated by RMD. However, these few genes had low -fold changes in general among both mock-infected and latency model comparisons, suggesting a possibility of enrichment in false positives in this gene subgroup.
We next sought to identify pathways that were affected by the HDACi treatment. Functional analysis of individual microarray expression (FAIME) (25) with RNA-Seq counts per million mapped reads as input and the KEGG pathway database were used for this analysis. Because there was a strong signal of differential expression in our data, many pathways were found to be affected by HDACi treatment. We have also compared these data with pathways perturbed in a different model of latency by ␣CD3/␣CD28 antibodies that induce T cell receptor (TCR) signaling (26), identifying multiple common pathways ( Fig. 4 and Table S2). Sixty-five pathways (25%) overlapped between the three treatments ( Fig. 4). Interestingly, subsets of the pathways in the overlap were modulated in opposite directions, as evidenced by the average difference in their FAIME scores (Table S2). This included metabolic pathways (e.g. genes in inositol phosphate metabolism pathway had higher expression following HDACi treatment and lower expression following TCR stimulation). Forty-three pathways (16.5%) overlapped between the two HDACis, 42 (16.1%) between RMD and TCR and 15 (5.7%) between SAHA and TCR (Fig. 4). Among these, identification of the estrogen signaling pathway specifically for RMD treatment and TCR stimulation, but not SAHA treatment, was a particularly interesting finding. The ESR1 gene encoding estrogen receptor 1 was down-regulated by RMD 24.3-fold, by TCR 25.6-fold, and much less strongly by SAHA (8-fold). ESR1 is required to maintain HIV provirus in a latent state (27); thus, its down-regulation by these treatments probably represents one of the secondary mechanisms of action by which latency is reversed. Interestingly, among genes modulated by RMD and not by SAHA, 17 were found to be members of the estrogen signaling pathway, consistent with the possibility that the specific effect on this pathway contributes to greater potency of RMD, compared with SAHA. EdgeR was used to identify differentially expressed genes by fitting general linear models and using the general linear model likelihood ratio test. Significance values were adjusted for multiple testing using the Benjamini and Hochberg method to control for FDR. Differentially expressed genes were filtered based on FDR Ͻ 0.05 and absolute -fold change Ͼ 1.5. Venn diagrams were constructed using the online tool Venny. Red, up-regulated genes; blue, down-regulated genes; darker shades represent overlapping sets of genes modulated by both SAHA and RMD.

Secondary mechanisms of action of HDACi and HIV latency A subset of genes modulated by SAHA and RMD has been implicated in regulation of HIV latency
We took advantage of the National Center for Biotechnology Information (NCBI) HIV-1 Human Interaction Database (28) to determine which of the thousands of genes modulated by HDACi were previously implicated in regulation of HIV latency. To identify genes with potential functions in the regulation of HIV transcription, genes from the NCBI database were searched against a list of keywords. The selected genes were then compared against the lists of genes modulated by SAHA and RMD (Fig. 5). The resulting gene lists were further curated manually by reviewing the Description column from the NCBI database to identify those genes modulated by HDACi that were the best candidates for regulation of HIV transcription. Among this reduced set of genes, 27 were up-regulated in common and 29 down-regulated in common by SAHA and RMD, whereas an additional 63 were modulated by RMD only (31 upand 32 down-regulated) and another 8 by SAHA only (6 up-and 2 down-regulated) ( Table 1 and Table S3). Among the genes up-regulated in common by both HDACis, 17 were associated with activation and 10 with repression of HIV transcriptional function. As noted above, RMD induced greater responses overall than did SAHA, with higher up-regulation in 15 of the 17 proposed "activators" and 6 of the 10 "repressors." Among the genes down-regulated by both HDACis, 15 were predicted to be repressors and 15 were predicted to be activators of HIV transcription (lysine acetyltransferase 5, KAT5, was represented in both groups due to its varied cell type-specific functions (29)). Again, RMD treatment produced greater responses than did SAHA, with increased down-regulation in 13 of the 15 repressors and 14 of the 15 activators. Among the genes that were modulated by RMD, but not SAHA, 18 repression-associated genes were up-regulated, and 15 activation-associated genes were down-regulated. TATA box-binding proteinassociated factor 10 (TAF10) and RB transcriptional corepressor 1 (RB1) were represented in both groups because both functions have been proposed in the literature (Table 1 and Table  S3). Altogether, these results indicate that RMD, like SAHA, exhibits multiple effects on host cell genes that would be predicted to either stimulate or inhibit HIV reactivation.

Genes predicted to stimulate or inhibit HIV reactivation from latency exhibit similar response kinetics to HDACi
In the RNA-Seq experiments, a single 24-h time point was analyzed, based on a previous protocol used to assess the activity of various LRAs (latency reversing agents) across different models of HIV latency and in study subject cells ex vivo (21). Because both HDACis showed multiple effects on gene expression associated with predicted inhibition of HIV reactivation, we took further steps to better delineate the timing of these potential inhibitory effects in relation to expression of those genes that are associated with stimulation of HIV reactivation. For this purpose, five genes that were modulated in common by SAHA and RMD were selected, based on their predicted roles in the regulation of early stages of HIV transcription (initiation and elongation), taken from their descriptions in the NCBI database. Genes with expected stimulatory effects on HIV reactivation included PR/SET domain 1 (PRDM1, also known as Blimp-1), a down-regulated repressor (30), and mediator complex subunit 26 (MED26), an up-regulated activator (31). Genes with expected inhibitory effects were represented by SWI/SNFrelated, matrix-associated, actin-dependent regulator of chromatin, subfamily B, member 1 (SMARCB1, also known as INI-1), a down-regulated activator (32), and poly(ADP-ribose) polymerase 1 (PARP1), an up-regulated repressor (33). Additionally, KAT5 (also known as Tip60) was selected, because its protein product has been reported to have cell-specific functions, as an activator in HeLa cells and a repressor in Jurkat cells (29).
First, we performed experiments to validate the observed changes in expression of the five selected genes (PRDM1, SMARCB1, KAT5, MED26, and PARP1) using a method independent from RNA-Seq. ddPCR was chosen for this purpose, as our laboratory has extensive experience with this technology in prior studies to validate gene expression (20,34). The expression of a housekeeping gene, ribosomal protein L27 (RPL27), which is not modulated by HDACi, was selected for normalization. The generated ddPCR data agreed well with the RNA-Seq data, with respect to direction and magnitude of the effects of SAHA and RMD on expression of the five selected genes. Differences in expression, between the HDACi-treated and DMSO-treated samples, were significant for both mock-infected cells and cells from the in vitro model of HIV latency (Fig. 6).
Because uninfected cells and cells from our HIV latency model responded to HDACi in the same manner for the five selected genes (Table S1), the kinetics experiments were performed using T cells from HIV seronegative, healthy volunteer blood donors, rather than setting up the in vitro model of HIV latency (n ϭ 6 and different from participants in the RNA-Seq study). CD4 ϩ T cells were treated with SAHA, RMD, or the solvent DMSO and examined after 6, 12, 24, and 36 h. ddPCR was used to measure expression of the five selected genes: PRDM1, SMARCB1, KAT5, MED26, and PARP1. After 24 h of treatment, there was very good agreement in the direction and magnitude of responses to SAHA and RMD, as measured by ddPCR, between the original sequenced data set and this new  26)) was conducted using FAIME. FAIME scores were obtained for each detected gene in each sample and compared using t tests. Significance values were corrected for multiple testing using the Benjamini and Hochberg method to control for FDR. Pathways with FDR Ͻ 0.01 were considered significant. The Venn diagram was constructed using the online tool Venny.

Secondary mechanisms of action of HDACi and HIV latency
data set. SAHA induced an early response, which achieved statistical significance (p Ͻ 0.001) at 6 h for all five selected genes. In contrast, the response to RMD was delayed; differences between RMD-and DMSO-treated samples achieved significance at 12 h and later (Fig. 7). Despite differences in the kinetics of response to SAHA and RMD, the time to peak of change of expression was similar both for up-and down-regulated genes, and for gene subsets, predicted to have stimulatory or inhibitory effects on HIV reactivation.

SAHA and RMD modulate expression of inhibitory genes in the same direction in the major maturation phenotypes of CD4 ؉ T cells
HDACi modulation of selected genes that are predicted to be inhibitory for HIV reactivation can be utilized to explore such gene pathways, as potential targets for improvement of shock and kill treatment strategies with HDACi. Of the five genes we selected, two were modulated in a manner that should inhibit HIV reactivation. As evidenced by the annotations from the NCBI database in the literature, the protein encoded by SMARCB1, a gene down-regulated by both SAHA and RMD, interacts with Tat and mediates Tat-dependent transactivation of the HIV promoter (32). The protein product of PARP1, a gene up-regulated by both SAHA and RMD, has a high affinity for the transactivation response element (TAR) of HIV RNA and binds to the TAR loop region to displace Tat or complexes of Tat with positive transcription elongation factor (P-TEFb) (33). Because KAT5 (down-regulated by HDACi) has demonstrated opposite effects in regulation of HIV transcription that depended on the targeted T cell line (29), we chose to focus on SMARCB1 and PARP1 for further evaluation.
We sought to determine whether all CD4 ϩ T cells that bear HIV proviruses are targeted uniformly by HDACi treatment. Before attempting to reverse any inhibitory HDACi effects, it is important to verify that the genes of interest are affected by HDACi treatment in the same manner in all of the relevant major maturation stages of CD4 ϩ T cells. For these experiments, we recruited eight new blood donors, who were differ-  (2), which were then used to annotate the differentially expressed genes (3). The final list was manually curated by evaluating the gene interaction descriptions in the NCBI database and relevant PubMed references, when necessary (4).

Table 1 Genes functioning as HIV activators or repressors that are modulated by HDACi
Function as HIV activator or repressor was determined from the NCBI HIV Human Interaction Database using the pipeline shown in Fig. 4. Listed are genes with absolute -fold change of Ն1.5 between HDACi and solvent DMSO treatment; all differences were significant (FDR Ͻ 0.05). TAF10, RB1, and KAT5 are listed as both activators and repressors because both functions were proposed in the literature (see Table S3 for details).

Secondary mechanisms of action of HDACi and HIV latency
ent from those who participated previously in the RNA-Seq and kinetics studies. Isolated CD4 ϩ T cells were sorted into the major maturation subsets that represent central memory, effector memory, and naive cells and were treated with SAHA (n ϭ 6) or RMD (n ϭ 5) for 24 h (three of eight donors had sufficient cell numbers to treat with both SAHA and RMD). Expression of SMARCB1 and PARP1 was mainly perturbed in the same direction (either up or down) in all three cell subsets (Fig. 8). One striking exception was a very weak and inconsistent up-regulation of PARP1 by RMD in the effector memory cell subset. In this case, the difference in up-regulation of PARP1 between effector memory and naive cells was statistically significant (p ϭ 0.0014; Fig. 8). There were no other significant differences in SMARCB1 and PARP1 modulation by the HDACi in the CD4 ϩ T cells of the major maturation phenotypes.

Discussion
Studies of secondary mechanisms by which LRAs may act to further enhance or, by contrast, inhibit HIV reactivation are warranted to develop alternative and more specific strategies for reactivation of latent HIV as well as to improve existing treatments by counteracting potential inhibitory effects. Analysis of the total transcriptome, following treatment with known LRAs, facilitates identification of modulated genes. In the present study, gene expression profiling was performed using RNA-Seq, a sensitive detection method for which analysis pipelines are well-developed. Genes modulated by two HDACis, SAHA and RMD, that have different potencies and HDAC specificities were evaluated.
In general, RMD is a more potent HDACi than SAHA. RMD functions at nanomolar concentrations to inhibit HDACs and Figure 6. Validation of gene expression modulation by HDACis using ddPCR. Same samples that were sequenced were subjected to ddPCR with assays to detect selected host genes. Data for both RNA-Seq and ddPCR are shown side-by-side for comparison. -Fold changes for RNA-Seq data were obtained by dividing read counts per million in HDACi-treated by read counts per million in DMSO-treated samples; -fold changes for ddPCR were obtained by dividing the corresponding number of RNA molecules detected in the ddPCR, normalized to housekeeping gene RPL27. Significance for ddPCR results was determined using a paired t test (*, p Ͻ 0.05; **, p Ͻ 0.01; ***, p Ͻ 0.001) for log 2 -transformed data, and all RNA-Seq results were significant (EdgeR's FDR Ͻ 0.05). Data are presented as mean -fold changes of gene expression following 24-h treatment with SAHA or RMD. Error bars, S.D.

Secondary mechanisms of action of HDACi and HIV latency
reactivate HIV from latency, whereas SAHA needs to reach hundreds of nanomolar to micromolar concentrations to show activity. RMD is predominantly active against HDACs of class I, HDAC1, -2, and -3, and significantly less active against HDACs from class II. SAHA, on the other hand, appears to be most active against HDAC6 of class II, although it also has activity against all class I HDACs (7). A recent study demonstrated that pan-HDACis, such as SAHA and RMD, induced expression of GFP reporter from HIV provirus in a greater number of cells, compared with class I selective HDACi (e.g. entinostat); in contrast, class I selective HDACis induced significantly more protein expression per cell. Moreover, when added to cells together with entinostat, SAHA and RMD reduced GFP reporter expression (19). A specific inhibitory mechanism was proposed for the case of SAHA, which may interfere with HIV expression via inhibition of HDAC6, thus causing changes in the acetylation state of its nonhistone targets, such as chaperone protein HSP90 (19), which is required for reversal of HIV latency (35)(36)(37). This mechanism of inhibition, however, would not be an expected mechanism for RMD, because it does not inhibit HDAC6. Nonetheless, RMD may disrupt the function of HSP90 chaperone complexes indirectly, via acetylation of HSP70, which caused a shift in association of the client proteins from HSP90 to HSP70, leading to their proteasomal degradation (38). Therefore, better understanding of the secondary mechanisms of action of LRAs, which may stimulate or inhibit HIV reactivation, is warranted for development of more potent shock and kill treatments.
Down-regulation of gene expression by both SAHA and RMD in our study is consistent with the presence of multiple secondary effects of both HDACis (i.e. gene modulation by other host factors induced as the result of chromatin relaxation

Secondary mechanisms of action of HDACi and HIV latency
following HDACi treatment). Interestingly, both up-and down-regulated genes had the same kinetics in our study (Fig.  7), consistent with a possibility that at least a subset of genes might be up-regulated not by initial chromatin relaxation but rather due to secondary effects induced by some of the initially up-regulated host factors. In our present study, these secondary effects of RMD were found to be much stronger than those of SAHA (greater magnitude of gene modulation by RMD as well as increased numbers of genes uniquely modulated by RMD only) (Table S1). Because RMD has higher overall potency to induce HDAC inhibition compared with SAHA, it is possible that its primary action affects more genes, and to a greater extent, which results in stronger and more diverse secondary effects.
Gene expression data from our study demonstrate that secondary mechanisms of action of both SAHA and RMD included modulation of expression of a number of genes implicated in regulation of HIV latency (Table 1 and Table S3). Importantly, effects on host genes that are predicted to be stimulatory for HIV reactivation occurred at the same time as the inhibitory effects for all five selected genes and both HDACis tested (Fig.  7), suggesting that complex interactions between HIV regulatory molecules likely exist to dictate the HIV reactivation outcome. An ability of RMD to reactivate HIV more potently than SAHA is likely due to the higher potency of RMD as an HDACi. However, specific effects of RMD on host genes and cellular pathways, such as ESR1 receptor and estrogen signaling, likely contribute to its greater potency as an LRA, compared with SAHA. Reactivation by both HDACi is likely hindered by the secondary mechanisms of action that may interfere with completion of transcription or induction of HIV protein expression. The net achievable level of induction of HIV gene expression likely depends on multiple host factors. Among samples analyzed in the present study, there was no relationship between percentages of cells bearing latent provirus and HIV reactivation (Fig. 2C). It is plausible that genomic locations of the proviruses, together with the host response to HDACi, contribute to the magnitude of response of HIV. The small sample size in our study did not allow us to perform correlation analyses of gene up-or down-regulation with HIV response, which warrants more detailed exploration in the future.
The HIV transcriptional profile in HIV-infected participants suggested that the extent of HIV transcriptional initiation is much greater than recognized previously and that the latent state is maintained by blocks to transcriptional elongation, completion, and splicing (39). Different LRAs exert differential effects on blocks of HIV transcription (39). In particular, treatment with RMD resulted in an increase of total and elongated transcripts but had less or no effect on polyadenylated and spliced transcripts (39). Both HDACis in our study also had the least effect on multiply spliced transcripts, despite the variability of responses among the participants (Fig. 2B). Combined with gene expression data, these results suggest that both SAHA and RMD likely have secondary effects on gene expression that may interfere with completion of HIV transcription. Whereas the present study focused on genes that mainly regulate HIV transcription, it is likely that host factors that positively or negatively regulate post-transcriptional stages of the HIV replication cycle are also modulated by HDACi (40).
The present study has begun to explore potential host targets that might be used to improve HIV reactivation in combination with HDACi to boost their ability to induce HIV gene expression by removing additional blocks to HIV transcription via reversal of the inhibitory secondary effects. Response to HDACi was highly variable, both in the present study (Fig. 2B) and in prior reports, in which cells from a subset of HIV-infected individuals were not responsive to SAHA (1), or a reduction of HIV RNA-positive cells was observed in a subset of cases following ex vivo treatment with RMD (41). One possibility is that CD4 ϩ T cells of various maturation states may respond to HDACi in a different manner (e.g. due to different availability of certain transcription factors (42) that may be modulated by HDACi). Thus, differences in responses to HDACi in various study participants may be explained by T cell composition. Therefore, before exploring a particular host gene as a potential therapeutic target to enhance HDACi-induced HIV reactivation, it is important to verify how the gene of interest responds to HDACi in cells of the relevant major maturation phenotypes that bear latent provirus. The results from our study demonstrate similar responses of SMARCB1 to HDACi in different maturation subsets (Fig. 8), suggesting that reversal of HDACi-induced SMARCB1 down-regulation would likely result in releasing the block to HIV transcription at the stage of RNA elongation in all cell types tested. The degree of PARP1 up-regulation by RMD depended on the maturation stage, with the least up-regulation induced in effector memory cells (Fig. 8). Thus, it is possible

Secondary mechanisms of action of HDACi and HIV latency
that interventions to counteract PARP1 up-regulation by RMD may be the least effective in effector memory cells.
One limitation of the present study was the reliance on NCBI database annotations to make predictions about the function of identified genes in HIV transcriptional regulation. Most of the studies referenced in the NCBI database were performed using cell lines or HIV reporter viruses, limiting transferability of these predictions to primary cells from HIV-infected individuals on ART. In light of some contradictory reports in the literature, such as in the cases of genes KAT5, TAF10, and RB1, validation of these findings in primary cells is very important. A recent study (43) has detailed the mechanism of action of KAT5 as an HIV transcriptional repressor and validated its role using a primary CD4 ϩ T cell model of latency and cells from HIVinfected individuals ex vivo. Inhibition of KAT5 activity with MG-149 enhanced HIV latency reversal induced with an LRA bromodomain inhibitor, JQ1 (43). Down-regulation of KAT5 expression by HDACi observed in the present study may represent one mechanism by which HDACi may further stimulate HIV reactivation following the initial chromatin relaxation. Another important example is PARP1. Whereas an older study cited in the NCBI database (33) proposed its function as an HIV transcriptional repressor via competition with Tat for binding TAR, a recent study (44) provided strong evidence that HIV transcriptional activation requires PARP1 activity. Different functional domains of PARP1 (DNA-binding and catalytic, respectively), have been implicated in these two studies. Because both studies were performed using cell lines only and neither used primary T cells, it would be particularly important to delineate the mechanism of action of PARP1 in relevant primary T cell model systems and to verify its role in latency reversal in cells from HIV-infected individuals. A model system that uses activated cells to establish latency and reporters to detect productively and latently infected cells (45) may be used to generate gene knockouts using a CRISPR/Cas9 system and to select for cells that contain both the gene knockout and HIV infection. One limitation of this approach, however, is low frequency of both the knockout and latently infected cells, resulting in insufficient cell numbers to perform subsequent treatments with HDACi and analysis. Because the CRISPR/Cas9 method requires T cell activation, it would not be applicable to models that use resting cells to establish latency or ex vivo cells from HIV-infected individuals. Another approach would be to use transfection methods; however, these are confounded by loss of cell viability, especially in the setting of extended primary cell culture necessary for latency establishment. We have currently begun to explore gene knockdown methods based on introduction of oligonucleotides called GapmeRs (Qiagen, Inc., Valencia, CA; previously Exiqon, Inc., Vedbaek, Denmark) into primary T cells or use of modified Accell siRNA (Dharmacon, Inc., Lafayette, CO), both of which are taken up by cells without transfection and might represent the method of choice for validating the role of selected genes in latency reversal in primary CD4 ϩ T cells.
An important observation of the present study is that in general, not many differences in response to HDACi due to the presence or exposure of cells to virus could be detected. Although the sample size was small, identification of 6 -11 genes as compared with thousands of genes modulated by each HDACi (Table S1) warrants the conclusion that the effect conferred by the presence of a provirus is negligible in latently infected cells. Therefore, future gene expression profiling studies of the effects of LRAs on primary CD4 ϩ T cells may be performed using uninfected cells to reduce time and cost of experiments. When the effect of the virus is of particular interest, large-scale single-cell strategies might be more appropriate.
Among genes differentially modulated by HDACi in the cells from the model of HIV latency and in mock-infected cells was CDKN1A/p21. Whereas this gene has been known to be modulated by HDACis and have a role in their anti-tumor properties (24,46), it was essentially not modulated at all in the in vitro model of HIV latency in the present study (Table S1). Because its expression was elevated in the cell line models of HIV latency, it was proposed as a latency biomarker (47). Its expression in elite controllers was higher, compared with progressors and HIV seronegative individuals; CDKN1A/p21 inhibited reverse transcription of the HIV genome and transcription from the integrated provirus via inhibition of cyclin-dependent kinase 9 (CDK9), a component of P-TEFb (48). In HIV-infected virologically suppressed individuals on ART, expression of CDKN1A/p21 negatively correlated with cell-associated HIV RNA and HIV transcriptional activity, defined as the ratio of HIV RNA to pol DNA (49). Altogether, the results from these studies are consistent with CDKN1A/p21 function to control HIV expression and maintain HIV in its latent state. Because HDACi did not up-regulate its expression in our latency model, CDKN1A/p21 is unlikely to interfere with induction of HIV expression by HDACi. However, evaluation of this effect in primary CD4 ϩ T cells from HIV-infected individuals would be needed to validate this finding.
In summary, our study demonstrates that SAHA and RMD modulate multiple host genes that have been implicated as HIV transcriptional regulators. Some of these genes may be explored as additional host targets for improving the outcomes of shock and kill strategies. The results from the present study suggest that complex interactions likely exist between HIV regulatory molecules to contribute to differential HDACi potencies and the outcome of HIV reactivation.

Study participants
Primary CD4 ϩ T cells from HIV seronegative volunteer blood donors were used for this study. The protocol was approved by the Institutional Review Boards of the University of California San Diego and the Veterans Affairs San Diego Healthcare System and abides by the Declaration of Helsinki principles. All volunteers gave written informed consent to participate in the study.

Isolation of primary CD4 ؉ T cells
CD4 ϩ T cells were isolated from whole blood using negative selection (StemCell Technologies, Inc., Vancouver, Canada). All CD4 ϩ T cell samples had Ͼ95% purity and Ͻ10% expression of activation marker HLA-DR, as assessed using the Accuri C6 flow cytometer (BD Biosciences). To isolate T cell maturation subsets, CD4 ϩ T cells were stained with ␣CD45RA antibody, con-

Primary CD4 ؉ T cell model of HIV latency
The primary CD4 ϩ T cell model of HIV latency (the Spina model) was set up as described previously (21), with minor modifications. Briefly, the experimental design (Fig. S1) involved setting aside a portion of the freshly isolated CD4 ϩ T cells ("bystander") and maintaining these cells in culture without infection or stimulation for 5 days. Another portion of CD4 ϩ T cells was stained with the viable dye e-Fluor 670, held overnight, and then infected with full-length HIV NL4-3 (4 -6 h) and stimulated on ␣CD3 ϩ ␣CD28 -coated plates. After 4 days of culture, the fully activated, infected, stained cells were mixed in with bystander cells at a ratio of 1:4 in the presence of the cytokines IL-2 and IL-15, to allow cell-to-cell virus transmission and establishment of latent infection in the bystander cells. At day 7 of culture following initial infection and activation, e-Fluor negative bystander cells were recovered by cell sorting (MoFlo XDP cell sorter, Beckman Coulter). On day 10, cells were treated with HDACi or the solvent DMSO in the presence of 0.5 M raltegravir (National Institutes of Health AIDS Reagent Program, catalog no. 11680) to prevent the spread of infection. A mock-infected cell culture condition was run in parallel with the set-up of the latency model (see legend to Fig. S1 for details) to determine whether exposure of cells to virus affects gene expression alterations induced by HDACi.

Treatment with HDACi
Aliquots of SAHA and RMD solubilized in DMSO were provided by Merck Research Laboratories, Inc. (Kenilworth, NJ) and Gilead Sciences, Inc. (Foster City, CA), respectively. Primary CD4 ϩ T cells or cells from the in vitro model of HIV latency were treated with SAHA (1 M), RMD (15 nM), or their solvent DMSO for 24 h. At harvest, cells in each sample were counted with a hemocytometer using trypan dye exclusion. For the follow-up kinetic studies, RMD was purchased from Selleckchem, Inc. (Houston, TX), and freshly resuspended aliquots were used for treatment. In the kinetics time course, cells were treated with HDACi and examined after 6, 12, 24, and 36 h. For studies of the effect of SAHA and RMD in CD4 ϩ T cells of different maturation states, cells were treated for 24 h. At collection, cells were counted using the Accuri C6 flow cytometer (BD Biosciences) and bead standards (Thermo Fisher Scientific).

Immunoblotting
CD4 ϩ T cells were treated with DMSO, SAHA, or RMD for 24 h. Cell pellets were lysed using radioimmune precipitation buffer containing protease inhibitors, and the concentration of protein was quantified using the Pierce BCA protein assay kit (Thermo Fisher Scientific). Ten micrograms of protein lysate was resolved on 4 -20% gradient SDS-PAGE and transferred to a polyvinylidene difluoride membrane (Bio-Rad) using a semidry transfer method. The membrane was probed with rabbit anti-acetylated histone H3 (MilliporeSigma, Inc. Burlington, MA) and goat anti-lamin B1 (Bio-Rad) primary antibodies and horseradish peroxidase-conjugated goat anti-rabbit and donkey anti-goat (Bio-Rad) secondary antibodies, respectively. Proteins were detected using the ChemiDoc MP Imaging System (Bio-Rad).

Flow cytometry analysis of T cell activation following treatment with HDACi
Cell activation was assessed in both the model of HIV latency and the mock-infected cells following treatment with HDACi. An aliquot of cells was stained with fluorophore-conjugated antibodies CD69-FITC clone L78, CD25-PE-Cy7 clone 2A3, HLA-DR-APC-Cy7 clone L243, and CD279-PE clone MIH4 and custom-made CD4-APC (BD Biosciences). The data were acquired using a FACSCanto II instrument (BD Biosciences) and analyzed using FlowJo version 10 software (FlowJo, LLC, Ashland, OR).

Integrant DNA assay
The integrant DNA assay was adapted from a previously published protocol (50). Total DNA was isolated using a Qiagen QIAamp DNA Micro Kit from latently infected cell pellets (0.5-1 million cells). The high-molecular weight fraction of genomic DNA (ϳ23 kb) was recovered from 0.5% low melt agarose gels using the Zymoclean Large Fragment DNA Recovery Kit (Zymo Research Corp., Irvine, CA). The DNA was then quantified using ddPCR for HIV gag, 2LTR circles, and RPP30 host cell genomic DNA (51). Total cell number was determined by dividing the number of detected RPP30 copies by 2. Integrated HIV DNA copy numbers were corrected by subtracting the 2LTR circles and normalized to the amount of cellular DNA input.

RNA-Seq data generation
ERCC spike-in RNA (Thermo Fisher Scientific) was added to cell lysates prior to RNA isolation (10 l of 1:800 dilution per million cells) to provide a means of normalization of resulting data for technical variations along the entire experiment. RNA was extracted using the Qiagen RNeasy kit. RNA integrity was assessed using a Bioanalyzer 2100 (Agilent Technologies) and deemed of good quality for conducting the RNA-Seq experiment; mean Ϯ S.D. RNA integrity numbers were on average 9.1 Ϯ 0.26. RNA-Seq was performed at Expression Analysis, Inc. (Morrisville, NC). RNA-Seq libraries were prepared using 100 ng of total RNA as starting material and the TruSeq Stranded Total RNA Library Prep kit (Illumina, San Diego, CA), after depletion of cytoplasmic and mitochondrial rRNA. All libraries were sequenced to a depth of 100 million reads using the Illumina HiSeq2500 sequencer to generate 50-bp pairedend library reads. Expression Analysis performed RNA-Seq data quality control, including filtering out low-quality reads, separating multiplexed samples by barcode, and removing 3Ј-adapter sequences, and provided resulting FASTQ files. FASTQ files are available through the Gene Expression Omnibus (GEO) accession number GSE114883.

RNA-Seq data analysis
To analyze expression of HIV, FASTQ files for each sample were mapped to the HIV NL4-3 genome (GenBank TM accession number AF324493) and counted using Bowtie (52). To analyze expression of host genes, reads were mapped to the GRCh38 human genome using TopHat (53) and counted against the human GENCODE (version 21) (54) using HT-Seq (55). All reads were mapped to the 92 ERCCs using Bowtie (52) and counted against individual ERCCs using HT-Seq. RUVSeq (56) was used to remove unwanted technical variation by using ERCC spike-ins. The EdgeR tool (57) in the R computing environment was used for differential gene expression analysis. EdgeR uses empirical Bayes estimation and exact tests based on the negative binomial distribution of RNA-Seq read data, followed by FDR correction using the Benjamini-Hochberg method (58). Genes were considered significantly modulated by HDACi when FDR-corrected p values were Ͻ0.05. To reduce the number of false positives, genes were also filtered by -fold change (absolute -fold change Ͼ1.5). For details on mapping, gene filtering, normalization, and linear modeling, see the supporting Methods. FAIME (25) was used to identify pathways affected by SAHA and RMD treatment and compare them with pathways modulated by ␣CD3/␣CD28 antibodies that induce TCR signaling. The latter data set was published previously (26) and is available from GEO under accession number GSE81810. Although ERCC spike-ins were also utilized in this data set, we have not used them for normalization in this analysis due to large shifts in ERCC spike-ins between groups associated with T cell activation. FAIME scores were generated using counts per million mapped reads as input data and compared between treatment and control conditions with paired t tests and FDR correction using the Benjamini-Hochberg method (58). FDR Ͻ 0.01 was considered a significant difference.
The NCBI HIV-1 Human Interaction Database (28) was then searched for genes that have been implicated in controlling HIV latency (see supporting Methods for details). EdgeR output was used to extract the expression information of genes of interest identified from the NCBI database. These resulting lists were manually curated to verify relevance to HIV latency, using the Description column of the NCBI database as well as available PubMed references. ddPCR ddPCR was performed as described previously (20). TaqMan assays for human genes were purchased from Applied Biosystems (now Thermo Fisher Scientific). The housekeeping gene, RPL27, was used to normalize expression of host genes. HIV RNA was measured using assays to detect unspliced (gag), polyadenylated (poly(A)), singly spliced (env), and multiply spliced transcripts. These assays were purchased from Integrated DNA Technologies, Inc. (Coralville, IA). HIV expression was normalized to cell number present in each sample based on live cell counts or measured copies of spike RNA ERCC-00051 added based on cell counts (Thermo Fisher Scientific). For assay IDs and sequences, see the supporting methods.

Statistical analyses
Flow cytometry data (percentages of cells bearing various activation markers) were compared across infection status and treatment conditions using ␤ regression modeling available in the betareg library (59) in the Bioconductor tool repository. To analyze ddPCR data, the -fold changes of expression following treatment with HDACi were obtained by dividing the RNA molecules per 1000 copies of RPL27 value in the HDACitreated sample by the DMSO-treated sample and log 2 -transformed. t tests were used to determine whether these log 2transformed -fold changes were different from 0 in ddPCR validation experiments. Time course data and comparison of gene modulation by HDACi among T cell subsets were analyzed using log 2 -transformed -fold changes for each treatment relative to DMSO at each time point and linear modeling (lm function) in R. All graphs were constructed using GraphPad Prism software (GraphPad Software, La Jolla, CA). Venn diagrams were constructed using the online tool Venny (J. C. Oliveros) and color-coded using Adobe Photoshop CS6.