S-adenosyl-methionine (SAM) alters the transcriptome and methylome and specifically blocks growth and invasiveness of liver cancer cells

S-adenosyl methionine (SAM) is a ubiquitous methyl donor that was reported to have chemo- protective activity against liver cancer, however the molecular footprint of SAM is unknown. We show here that SAM selectively inhibits growth, transformation and invasiveness of hepatocellular carcinoma cell lines but not normal primary liver cells. Analysis of the transcriptome of SAM treated and untreated liver cancer cell lines HepG2 and SKhep1 and primary liver cells reveals pathways involved in cancer and metastasis that are upregulated in cancer cells and are downregulated by SAM. Analysis of the methylome using bisulfite mapping of captured promoters and enhancers reveals that SAM hyper-methylates and downregulates genes in pathways of growth and metastasis that are upregulated in liver cancer cells. Depletion of two SAM downregulated genes STMN1 and TAF15 reduces cellular transformation and invasiveness, providing evidence that SAM targets are genes important for cancer growth and invasiveness. Taken together these data provide a molecular rationale for SAM as an anticancer agent.


SUPPLEMENTARY MATERIALS Cell culture and drug treatments
HepG2 and SkHep1 cells were maintained in MEM medium (Gibco), supplemented with 2 mmol/L glutamine (Sigma-Aldrich), 10% fetal bovine serum (FBS; Gibco), 1 U/mL penicillin, and 1 μg/mL streptomycin (Gibco). NorHep cells were maintained in human hepatocyte cell culture complete medium (Celprogen). Cells were grown in a humidified atmosphere of 5% carbon dioxide at 37°C. Following 3 to 5 minute incubation with Trypan blue after trypsinization, the viable cells were counted under the microscope.
The ability of cells treated with SAM to invade through extra cellular matrix was evaluated by the Cell Invasion Assay Kit (Chemicon Int.). The kit utilizes a reconstituted basement membrane matrix of proteins derived from Engelbreth-Holm-Swarm (EHS) mouse tumor. Briefly, 50,000 cells resuspended in serum-free media were added to the inserts dipped in the lower chamber containing complete media. Following 24h incubation at 37°C, invasive cells were stained and counted under the microscope. Additionally, 50,000 viable cells (as determined by Trypan blue) resuspended in complete media were added to a six-well plate and were counted following 24h incubation period concurrently with measuring invasiveness to measure cell viability.
To determine anchorage-independent growth on soft agar, a measure of transformation in vitro [1], 3,000 viable cells treated for 5 days with SAM were seeded into soft agar and plated in triplicate in a 6-well dishes containing 4mL of complete medium with 0.33% BD Bacto™ agar solution at 37°C as previously described [2]. The total number of colonies (>10 cells/colony) that formed on soft agar was counted under the microscope after three weeks of plating.

DNA/RNA extraction, quantitative real-time PCR and western blot
DNA and RNA was extracted using AllPrep DNA/ RNA/miRNA Universal Kit (Qiagen) according to the manufacturer's protocol. One microgram of total RNA served as a template for cDNA synthesis using 20U of AMV reverse transcriptase (Roche Diagnostics), as recommended by the manufacturer. The quantitative real-time PCR (QPCR) reaction was carried out in Light Cycler 480 machine (Roche) using forward and reverse primers listed in Supplementary Table 1. Quantification was performed using Roche LightCycler 480 software second derivative method.

shRNA inhibition
For STMN1 and TAF15 depletion we used the lentivirus human pGIPZ shRNA plasmids and control pGIPZ-scrambled shRNA (Open Biosystems) (Supplementasy Table S2). Lentiviruses were assembled using the following three vectors: GFP expression pGIPZ transfer vector-includes the insert (Open Biosystems); pMD2.G (VSV-G envelope expressing plasmid); PAX (packaging plasmid). The day before transfection, 10 6 HEK293T cells were plated in a 10 cm dish (20-30% confluence). Next day, 5 μg of each vector were transfected using FuGene HD transfection reagent (Roche) according to the manufacturer's protocol. Cells were incubated for 48 h followed by collection of the medium containing the virus, filtered and used to infect the target cells. Selection with 1 mg/ml puromycin (Sigma) was started 48h post infection. Specific shRNAs were selected based on knockdown efficiency in the specific cell line; STMN1 was targeted with ShSTMN1#V3LHS_383505 and TAF15 was targeted with ShTAF15#V2LHS_172493 (Supplementary Table 1 for sequences).

RNA sequencing and data analysis
RNA (4 μg) from SAM and SAM buffer (control) treated cells was processed using the TruSeq RNA Sample Prep Kit (Illumina, San Diego, USA) following the manufacturer's protocol. Briefly, polyA mRNA was purified and fragmented. Then first and second strand cDNA synthesis was performed, ends were repaired and adenylated, adapters were ligated, and fragments were enriched with PCR amplification. The library was validated with an Agilent Technologies 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA). RNA-seq preparations were sequenced using Illumina HiSEQ2K platform (Illumina, San Diego, USA). A total of 50bp pair-end reads were generated. Duplicate samples were sequenced. Around 30 million reads were obtained per sample (Supplementary Table 3).
Fastqc was used to QC RNAseq data. Raw reads were trimmed sequentially for adapter contamination. The paired end reads were aligned to the human reference sequence (hg19, Feb. 2009) with TopHat 2.0.9 [4] with default setting, except for the option -g. Aligned bam files were assessed using cufflinks v.2.2.1 [4]to generate a transcriptome assembly and to estimate the expression level (FPKM) of all detected isoforms. FPKM was calculated as number of end-paired reads (a single fragment per end-paired reads) mapped to a gene divided by the number of all fragments mapped to the genome (in million) and the length of RNA (in KB). The tagwise dispersion were estimated and then used for logFC (log 2 fold change) estimating and testing. Differentially expressed genes (DEGs) were extracted by applying the threshold false discovery rate (FDR) of less than 0.05 to adjusted P values using Fisher's exact test. Furthermore, Reads counts were obtained from the mapping results by using HTSeq-count (1.0) [5]. Extreme low expressed genes < 2 count-per-million (CPM) were filtered using EdgeR package [6].

DNA capture bisulfite sequencing
DNA (5μg) from SAM and SAM buffer (control) treated cells was used for preparation of DNA capture bisulfite sequencing library. The DNA libraries were prepared using Illumina's paired-end sequencing DNA sample preparation kit according to the manufacturer's protocol. Target DNA fragments were captured using the human SeqCap Epi CpGiant Enrichment Kit (Nimbergen, USA) according to the manufacturer's recommendations. SeqCap Epi CpGiant (Nimblegen, Roche) interrogates more than 5.5 million CpG sites covering promoters and regulatory sequences in human genome. In order to capture bisulfite converted DNA, probes are designed to hybridize to both strands of fully methylated, partially methylated and fully unmethylated derivatives of the genomic target and then pooled together. DNA libraries were qualitychecked and quantified on Agilent 2100 Bioanalyzer. Sequencing was performed on the Illumina HiSEQ2K platform using a standard 50 cycle paired-end read sequencing protocol and Illumina's sequencing reagents according to the manufacturer's recommendations.

Bisulfite sequencing data analysis
The raw data was processed as recommended by Sequencing Solutions Technical Note from Roche for SeqCap Epi CpGiant bisulfite sequencing data analysis (https://sftp.rch.cm//diagnostics/ sequencing/literature/nimblegen/07292163001_NG_Seq Cap_TchNote_EvalEpiData.pdf). Briefly, the FASTQ reads were trimmed using Trimmomatic0.30 [7] and aligned with BSMAP2.74 [8] to the human reference sequence (hg19, Feb. 2009). After removal of duplicated reads and filtering for properly paired reads from BSMAP mapping results, methylation ratios were extracted using python script methratio.py from the BSMAP package. Methylation difference was calculated using methylKit [9] according to user guide (coverage count >5). The differentially methylated CpGs were extracted with a q-value <0.05 and delta methylation >15%. The differentially methylated sites were annotated with CHIPseeker package [10].

Statistical analysis
Statistical analysis of pyrosequencing, QPCR, invasion, soft agar and cell viability assays was performed using an unpaired t test with two tailed distribution. The results were considered statistically significant when P <0.05.
The Ingenuity Pathway Analysis (IPA) program (http://www.ingenuity.com/index.html) was used to compute enriched gene networks, functional categories, canonical pathways and upstream regulators. Heatmaps were created using GeneE (http://www.broadinstitute.org/ cancer/software/GENE-E/doc.html). GSEA (GSEAv2.2.3) [11] was used to examine the gene sets that are enriched with genes whose expression increases after SAM treatment that are significantly downregulated in cancer vs. normal cells and genes whose expression decreases after SAM treatment that are significantly upregulated in cancer vs. normal cells. GSEA was also used to examine the gene expression profiles of gene sets enriched with genes that are either hypomethylated or hypermethylated in response to SAM. For differential expression FDR of < 0.05 and countper-million (CPM) of >=2 were used as thresholds. For differential methylation the thresholds were q-value <0.05 and delta methylation >15% in promoter regions (-1000bp-1000bp). Significantly enriched gene sets after 1,000 permutations at FDR of <0.25 were reported.