SKM used for these studies was obtained from animals used in experiments approved by the Institutional Animal Care and Use Committee at both Tulane National Primate Research Center (TNPRC) in Covington, Louisiana and Louisiana State University Health Sciences Center (LSUHSC) in New Orleans, Louisiana, and adhered to National Institutes of Health guidelines for the care and use of experimental animals. The pathophysiological course of SIV infection has been previously described in published manuscripts [11–14]. A total of 28 4–6-year-old male macaques (Macaca mulatta) were studied in three experimental groups: SUC/SIV group (n = 9), CBA/SIV group (n = 11), and control group (n = 8). For the control group, skeletal muscle samples were obtained at necropsy from a group of SIV-negative, healthy control macaques and used as reference values for comparison of analyzed variables.

Experimental protocol

The gastric catheter placement for alcohol delivery, the alcohol delivery protocol, and the route of SIV infection have been described in detail elsewhere [11, 14, 32]. Animals were briefly exposed to daily intragastric administration of a mean of 2.5 g per kg body weight ethanol (30 % w/v water), beginning 3 months prior to SIV infection and continuing throughout the duration of study. This protocol of CBA administration provided an average of 15 % of the animals’ total daily caloric intake and produced blood alcohol concentrations of 50–60 mM. Control animals were infused with sucrose. Animals were provided monkey chow (Lab Fiber Plus Primate diet DT; PMI Nutrition International, St. Louis, MO) ad libitum and supplemented with fruits, vitamins, and Noyes treats (Research Diets, New Brunswick, NJ).

Three months after initiating CBA administration, animals were inoculated intravenously with 10,000 times the 50 % infective dose (ID 50 ) of SIV mac251 . SIV disease progression was monitored throughout the study period through clinical, biochemical, and immunological parameters (CD4/CD8 lymphocyte ratios) in addition to plasma viral kinetics (SIV gag RNA levels) as described elsewhere [12, 32]. SKM (gastrocnemius) samples were obtained at necropsy when animals met any one of the criteria for euthanasia based on the following: (1) loss of 25 % of body weight from maximum body weight since assignment to protocol, (2) major organ failure or medical conditions unresponsive to treatment, (3) surgical complications unresponsive to immediate intervention, or (4) complete anorexia for 4 days. SKM tissue samples were dissected, snap frozen, and stored at −80 °C until analyses. SKM samples from all animals in the three treatment groups were used for the gene microarray, microRNA microarray, methylation array and qPCR. SKM samples used for the analysis in this study were used in two previously published studies [13, 14].

mRNA Microarray analysis

To determine how CBA affects the expression of mRNA at end-stage SIV infection, total RNA was isolated from SKM at necropsy from a total of 28 animals (SUC/SIV (n = 9), CBA/SIV (n = 11), and control (n = 8)) and microarray analysis was performed as previously described on all samples [13]. Briefly, the microarray hybridization was performed at the Stanley S. Scott Cancer Center’s Translational Genomics Core at LSUHSC in New Orleans, Louisiana. Total RNA was extracted using the RNeasy Mini Kit (Qiagen, Valencia, CA) according to the manufacturer's instructions. The RNA was hybridized to Illumina HumanWG6_v3 chips (San Diego, CA) following manufacturer’s instructions. For data analysis the samples were normalized using the cubic spline algorithm, assuming that the distribution of transcripts is similar [33]. Differential expression was determined by comparing treatment and control groups using the Illumina Custom algorithm that assumes that target signal intensity is normally distributed among replicates corresponding to some biological condition.

The fold change in gene expression of SUC/SIV and CBA/SIV was obtained by dividing the expression level of each gene by that of the control. The fold change of CBA/SIV/SUC/SIV was obtained by dividing the expression level of each gene between CBA/SIV and SUC/SIV. Comparisons between the gene expression levels of CBA/SIV with SUC/SIV animals reflect the impact of CBA on SIV-mediated changes in gene expression. Genes whose expression was altered in an alcohol-dependent manner by ≥1.5-fold were examined. The GEO accession number is GSE59111.

Gene Set Enrichment Analysis (GSEA) was run on the normalized, unfiltered microarray dataset as suggested in the tools implementation (http://software.broadinstitute.org/gsea/index.jsp) version 2.2.1) [34, 35] using the c5.all.v5.symbols. (Gene ontology), running 1000 permutations and excluding gene sets with fewer than 5 genes or more than 500. Using these parameters, 111 gene sets were selected for the analysis. The GSEA statistics are detailed in http://www.broadinstitute.org/gsea/doc/GSEAUserGuideFrame.html [36].

MicroRNA microarray analysis

To determine the impact of CBA on miRNA expression at end-stage SIV infection, small RNA (<30 bp) were purified. The microRNA microarray hybridization was performed at the microarray core at LSUHSC. Briefly, 500 μg of total RNA was biotin-labeled using the FlashTag Biotin HSR kit (Genisphere, Hatfield, PA) according to the manufacturer’s instructions. All samples showed expected labeling and the resulting targets were hybridized to miR 2.0 arrays (Affymetrix) containing 15,644 microRNA probe sets from the miRBASE v15 (http://microrna.sanger.ac.uk), which recognize microRNAs from a number of organisms including human and rhesus macaques. The arrays were washed and processed on a Fluidics Station 450 and scanned with a confocal laser scanner (GeneChip Scanner 3000, Affymetrix) according to the manufacturer's instructions. Data from the microarrays were analyzed with Affymetrix microRNA QC Tool according to the manufacturer’s instructions. Differential miR expression was performed using the One-Way ANOVA workflow with multiple test correction (FDR = 0.05) and selected from volcano plot using p-value <0.05 and exhibiting a ≥ 1.5-fold change cutoff.

All differentially regulated miRs were analyzed using miRSystem (http://mirsystem.cgm.ntu.edu.tw/) [37]. miRSystem is a database which integrates seven miRNA target gene prediction programs: DIANA, miRanda, miRBridge, PicTar, PITA, rna22, and TargetScan. The database also contains validated data from TarBase and miRecords. To balance the reliability of the predictions, three algorithm hits were used: 1) hypergeometric p-value determination; 2) empirical p-value is determined by ranking the enriched hypergeometric probability as compared with null baseline probabilities. The null baseline probability was established by randomly selecting a group of miRNAs, between 1 and 100, and using the default values in mirsystem to calculate the raw p-value for each pathway [37]; and 3) a weighted pathway-ranking method is also calculated from the expression ratios of the differentially regulated miRNAs to rank the enriched pathways. For each functional category, the ranking score was obtained by summation of the weight of its miRNA times its enrichment -log (p-value) from the predicted target genes \( \left\{\mathrm{Score}={\displaystyle \sum \mathrm{AmiRNAwi}\left[\hbox{-} \log 10\left(\mathrm{pi}\right)\right]}\right\} \). The target genes, pathway ranking, and functional annotation summaries are included in the results. The parameters set for analysis were: a) validated genes equal to or more than 3, b) observed to expected ratio (O/E) greater or equal to 2 and c) total genes in the pathway ≥ 25 and ≤ 500.

We also analyzed individual differentially expressed miRNAs using miRTarBase and TargetScan. A search for target genes using miRTarBase [38] determined that 4 of the downregulated miRNA target 12 independent genes; targeting is validated in the literature by at least two independent experimental methods [38]. For the remaining miRNAs that have no experimentally validated gene targets, we used the TargetScan database, which identifies genes statistically predicted to be targets, using Homo sapiens as the reference. This search identified 8 individual downregulated miRNAs predicted to target 18 independent genes based on their predicted efficacy of targeting (context score) [39, 40] or probability of conserved targeting (P CT ) [41] [context score ≥ 85 %; P CT ≥ 0.8] [42].

DNA methylation microarray analysis

To determine how CBA affects promoter methylation in SIV-infection, the Infinium HumanMethylation27 array was utilized. The Infinium HumanMethylation 27 examines more than 27,000 CpG islands in more than 14,000 genes’ promoters. SKM samples were bisulphite-converted with Zymo EZ DNA Methylation kit (Zymo Research, Irvine, CA, USA). GenomeStudio v2011.1 (Illumina, San Diego, CA, USA) with Methylation module 1.9.0 software was used in the methylation analysis. The Infinium platform assays covers 96 % of CpG islands with multiple sites in the island, the shores (within 2 kb from CpG islands), and the shelves (>2 kb from CpG islands). All the Illumina quality controls were acceptable, including sample-independent and dependent controls, staining controls, extension controls, target removal controls, hybridization controls, bisulphite conversion I and II controls, specificity controls, non-polymorphic controls and negative controls. Probes were considered to be differentially methylated if the resulting adjusted p-value was <0.05. The Benjamini-Hochberg method [43] was used to adjust the p-values and ensure that the false discovery rate was <0.05. The corresponding gene list was derived from the gene annotations associated with the probes. The GEO accession number is GSE75729. Functional enrichment analyses of genes with differential methylation of promoter regions were performed using DAVID Bioinformatics Resources. The gene ontology Biological Processes (BP_all) is represented for genes. The count: number of genes involved in the term; %: percentage of involved genes/total genes; p-value: modified fisher exact p-value, EASE Score (p ≤ 0.05); Benjamini: adjusted p-value using Benjamini-Hochberg procedure is presented in the results [44, 45].

The biological functions of individual differentially expressed genes in the mRNA microarray, validated target genes of differentially expressed miRNAs and genes altered in the DNA methylation array were also categorized by searching each individual gene in the GeneCards® Human Genome Database (http://www.genecards.org). This is an integrative database that provides descriptions of gene functions as extracted from multiple public databases, including Entrez Gene, UniProtKB, Tocris, Bioscience, and PharmGKB. Additional information was derived from searching each individual gene in the NCBI Gene Database (http://www.ncbi.nlm.nih.gov/gene/); the remaining genes with unclear biological function were analyzed using the DAVID Bioinformatics Database (https://david.ncifcrf.gov). This allowed for analysis of genes with specific functions relevant to skeletal muscle function. Those genes that did not fall within any obvious biological function after these analyses were classified as “Miscellaneous.” This analysis is included as supplementary tables.

q PCR (qPCR) for miR expression

To validate the deep sequencing data, the relative expression of 3 differentially expressed miRs (miR-34a, miR-10b, miR-20) was further determined by individual Taqman miR assays (Thermo Fisher Scientific). Approximately 200–250 ng of total RNA was reverse-transcribed using the stem loop primers provided in the predesigned kit and ~1.3 μl of cDNA was subjected to 40 cycles of PCR on the CFX96 Bio-Rad PCR cycler (Bio-Rad) using the following thermal cycling conditions: 50 °C for 2 min, 95 °C for 10 min followed by 40 repetitive cycles of 95 °C for 15 s and 60 °C for 1 min. As a normalization control for RNA loading, SNOU6 were amplified in duplicate wells on the same multi-well plate.

qPCR for target genes of differentially regulated miRNAs

Total RNA isolated for miR sequencing studies was used for gene expression analysis as well. cDNA was synthesized from 1000 ng of the resulting total RNA using the Quantitect Reverse Transcriptase Kit (Qiagen), in accordance with the manufacturer’s instructions. Primers were designed to span exon-exon junctions (IDT, Coralville, IA) and used at a concentration of 500 nmol. The final reactions were made to a total volume of 20 μl with Quantitect SyBr Green PCR kit (Qiagen). All reactions were carried out in duplicate on a CFX96 system (Bio-Rad Laboratories, Hercules, CA) for qPCR detection. qPCR data were analyzed using the comparative Ct (delta-delta-Ct, ΔΔCT) method. Target genes were compared with the endogenous control (ribosomal protein S13 (RPS13)) and CBA/SIV and SUC/SIV values were normalized to controls.

Availability of supporting data

The data sets supporting the results of this article are available in Gene Expression Omnibus (GEO). The GEO accession numbers for the microarray data are GSE59111 and GSE75729. The links to the data are http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE75729 and http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE59111.