Abstract
Background/Aim: Vascular anomalies encompass different vascular malformations [arteriovenous (AVM), lymphatic (LM), venous lymphatic (VLM), venous (VM)] and vascular tumors such as hemangiomas (HA). The pathogenesis of vascular anomalies is still poorly understood. Viral infection was speculated as a possible underlying cause. Materials and Methods: A total of 13 human vascular anomalies and three human skin control tissues were used for viral analysis. RNA derived from AVM (n=4) and normal skin control (n=3) tissues was evaluated by RNA sequencing. The Virome Capture Sequencing Platform for Vertebrate Viruses (VirCapSeq-VERT) was deployed on 10 tissues with vascular anomalies (2×AVM, 1×HA, 1×LM, 2×VLM, 4×VM). Results: RNA sequencing did not show any correlation of AVM with viral infection. By deploying VirCapSeq-VERT, no consistent viral association was seen in the tested tissues. Conclusion: The analysis does not point to the presence of an active viral infection in vascular anomalies. However, transient earlier viral infections, e.g. during pregnancy, cannot be excluded with this approach.
Vascular anomalies (VA) of the head and neck area are a highly heterogeneous group of rare vascular diseases (1-3). VA can have a major clinical impact since, depending on their size and anatomical location, patients can develop difficulties in breathing and swallowing. Furthermore, and particularly true for high-flow vascular malformations, such as arteriovenous malformations (AVM), patients can suffer from extensive and life-threatening bleeding. The pathogenesis of VA is still poorly understood, despite recent findings pointing to mutations in genes such as phosphatidylinositol-4,5-bisphosphate 3-kinase catalytic subunit alpha (PIK3CA) (4) and angiopoietin receptor (TIE2) (5) as a possible underlying cause. Many VA are already present at birth, indicating that the pathological impact leading to disturbance of the vascular architecture, as observed in VA, has occurred during pregnancy. For years, viral infections have been suspected as one possible cause for the development of VA (6). However, lack of advanced techniques allowing for comprehensive viral detection made it difficult to come to a clear conclusion. A previous study discussed the impact viruses can have on the vasculature of the body (7). Interestingly, it was observed that endothelial cells and neural structures in the vicinity of human papillomavirus (HPV)-positive cervical cancer appeared positive for viral proteins and genes (8). Similarly, our group previously noted a possible association of HPV infection with VA during immunohistochemical analysis (unpublished observations), which, however, were not confirmed in subsequent studies. This prompted us to search for methods allowing for comprehensive viral detection in tissues with VA, taking into consideration not only HPV but also other human pathogenic viruses. For many years, studies aiming to evaluate a possible association of human diseases and viral infection have suffered from the restriction of detection methods available. Most commercial detection methods were established only for the detection of single viruses or virus families, thereby not allowing recording of several different viruses in parallel. The development of high-throughput sequencing techniques now allows for a comprehensive evaluation of changes in the genome and transcriptome of cells (9).
In the present study, we evaluated the presence of viral transcripts in human VA using RNA sequencing (RNA-Seq) and the highly sensitive Virome Capture Sequencing Platform for Vertebrate Viruses (VirCapSeq-VERT) (10).
Materials and Methods
Tissues. All VA [arteriovenous (AVM), lymphatic (LM), venous lymphatic (VLM), venous (VM), and hemangiomas (HA)] and control tissues were obtained from patients treated at the Department of Otorhinolaryngology (University Hospital Marburg, Germany) according to the guidelines of the local Ethics Committee (reference number: 111/12; Medical School, Philipps University Marburg, Germany) and in cooperation with the local Biobank of the Medical School (Director: Professor P.I. Pfefferle). A total of 13 human VA and three human skin control tissues were used for viral analysis. RNA derived from AVM (n=4) and normal skin control (n=3) tissues was evaluated by RNA-Seq. VirCapSeq-VERT was deployed on 10 VA tissues [2×AVM, 1×LM, 2×VLM, 4×VM, and 1× HA]. Details regarding the tissue specimens used in the study are summarized in Table I.
Cottontail rabbit papillomavirus-associated VX2 tumor of rabbit served as a positive control. The animal experiment was approved by the Regional Board Giessen, Germany [V54-19c20-15(1) MR20/26 Nr.34/2011] according to the German Animal Protection Law. Auricular tumors were experimentally induced by inoculation of a VX2 tumor cell suspension as described previously (11). Auricular VX2 tumor samples were randomly taken from two adult male Iffa Credo New Zealand White rabbits (Charles River Wiga Deutschland GmbH, Sulzfeld, Germany).
Immunohistochemistry. Tissues were formalin-fixed and paraffin-embedded (FFPE). Immunohistochemistry was performed according to a standard laboratory protocol as reported elsewhere (12). In short, 3 μm-thick FFPE slices were generated using a microtome followed by deparaffinization (Rotihistol; Carl Roth, Karlsruhe, Germany), alcohol incubation and heat treatment in 0.01 M tri-sodium citrate dihydrate buffer (92-95°C, 20 min). Samples were blocked in goat normal serum and incubated (1:100) with anti-CD31 specific antibody (Dako, Glostrup, Denmark) followed by incubation with the anti-mouse Envision+ System - HRP Labelled Polymer (Dako). Antibody polymer complex binding was visualized after incubation with 3,3’-diaminobenzidine substrate (Dako). Mouse IgG (Dako) was used as a negative control and did not exhibit any significant background (not shown). Samples were counterstained with Mayer's Hemalaun solution (Merck, Darmstadt, Germany). Microscopic documentation of staining results was performed with a Leica DM IL LED microscope (Leica Microsystems, Wetzlar, Germany) using Leica Application Suite software (version 4.8.0).
RNA isolation and validation. For RNA-Seq, human AVM (n=4) and skin control (n=3) tissues as well as rabbit VX2 tumors (n=2) (Table I) were homogenized with Precellys 24 Homogeniser (Bertin GmbH, Frankfurt, Germany). The homogenate was used to extract whole cellular RNA (RNeasy® Mini Kit; Qiagen, Hilden, Germany). RNA quantity and quality were evaluated with the NanoDrop™ 2000 photometer (Thermo Fisher Scientific, Darmstadt, Germany) and the Experion™ Automated Electrophoresis System (Bio-Rad Laboratories GmbH, München, Germany) prior sending it to the European Molecular Biology Laboratory Genomics Core Facility (Heidelberg, Germany). Total RNA from normal adult rabbit skin was obtained from BioCat GmbH (BioCat GmbH, Heidelberg, Germany). For VirCapSeq-VERT, VA tissues (n=10; see Table I) were homogenized followed by nucleic acid extraction with the Qiagen All Prep Kit (Qiagen) and quality check.
RNA-Seq. Qubit RNA BR Assay kit (Q10211; Invitrogen™, Thermo Fisher Scientific) was used to assess the quality and quantity of total RNA. RNA integrity was evaluated with the 2100 Bioanalyzer System (Agilent Technologies, Santa Clara, CA, USA) using RNA 6000 Pico Reagents and Chips according to the manufacturer's instructions. An RNA Integrity number of 7 or higher was the criterion to proceed with library preparation. A total of 300 ng of validated total RNA in 50 μl water was deployed for library preparation (all samples were normalized to the same amount of RNA). Library preparation was performed with the TruSeq® Stranded mRNA Library Prep Kit (Illumina Inc., San Diego, CA, USA). In short, polyA mRNA was isolated from total cellular RNA using Oligo (dT) attached magnetic beads. Poly-A RNA was fragmented and cDNA first-strand synthesis (25°C for 10 min; 42°C for 15 min; 70°C for 15 min; 4°C, hold) was performed in the presence of hexamer primers followed by second-strand synthesis (Second Strand Marking Master Mix; Illumina Inc.) for 1 h at 16°C. DNA purification was carried out with SPRISelect beads (B23319; Beckman Coulter GmbH, Krefeld, Germany) in the same manner as described in the Illumina protocol for AMPure XP beads. Purified dsDNA was eluted from the beads in 15 μl of resuspension buffer. The A-tailing reaction was performed at 37°C for 30 min followed by 5 min heat inactivation at 70°C. Subsequently, ligation of adapters (containing a unique 7-bp sequence for indexing) was implemented for all 10 samples (each sample received a unique adapter). After ligation, the DNA was purified and size selected according to the instructions of Illumina using SPRISelect beads. A polymerase chain reaction step was included to enrich DNA fragments with adapter molecules on both sides (step 1: 98°C for 30 s; step 2: 98°C for 10 s; step 3: 60°C for 30 s; step 4: 72°C for 30 s; steps 2-4 ×15; step 5: 72°C, 5 min). After a final purification step (SPRISelect beads), the quantity and quality were assessed with Qubit ds DNA HS Assay Kit (Invitrogen™). The fragment profile as well as the presence of adapter dimers were analyzed with an Agilent 2100 Bioanalyzer System (Agilent DNA 1000 Reagents and Chips; Agilent Technologies). The barcoded libraries were pooled at equimolar concentrations followed by dilution to 4 nM prior sequencing. Five microliters of the diluted pool were combined with 5 μl of 0.2 N NaOH for denaturation (5 min at room temperature) followed by neutralization with 5 μl Tris-HCl (200 mM, pH 7.0) and addition of 985 μl Illumina HT1 buffer. The sample (20 pM) was further diluted in HT1 buffer to yield a final loading concentration of 1.8 pM. The sample was then run on an Illumina NextSeq500 instrument in the 75 + index read cycle mode. RNA-Seq was performed at the Genomics Core Facility at the European Molecular Biology Laboratory.
Tissues used for viral transcript analysis in this study.
For statistical analysis, the unpaired Student's t-test was used to evaluate differences between human skin and human AVM, respectively. Welch's correction was applied whenever the variances of the group means differed significantly.
VirCapSeq-VERT. Total cellular RNA was isolated from 10 VA tissues (Table I) as described above. The subsequent procedure was performed according to an established protocol at the Center for Infection and Immunity (Columbia University, New York City, NY, USA) (10). In short, total cellular RNA was reverse transcribed into cDNA with random priming (Superscript™ III; Invitrogen™) followed by second-strand synthesis (Kleenow, New England Biolabs, Frankfurt/Main, Germany). The double-stranded DNA was sheared to an average of 250 bp with a Covaris E210 focused ultrasonicator followed by fragment size evaluation (Bioanalyzer 2100; Agilent). Subsequently, double-stranded DNA fragments were subjected to DNA end repair, dA-tailing and indexed adapter ligation followed by PCR-amplification (Hyper Library Preparation kit, KAPA; KapaBiosystems, Wilmington, MA, USA). Amplicons were then hybridized with the biotinylated VirCapSeq-VERT capture probes. Probe-bound cDNAs were trapped by magnetic streptavidin-beads and washed (SeqCap EZ Reagents; Roche NimbleGen, Inc., Madison, WI, USA). After low-cycle PCR amplification, the libraries were sequenced on the HiSeq 2500 system (Illumina; 10 index codes per lane).
Bioinformatics. RNA-Seq data were provided for further analysis after quality control and mapping to the human reference genome hg19. Files for rabbit samples were converted back into fastq-files using SAMtools (13) and mapped against the Oryctolagus cuniculus reference genome (GenBank Assembly ID GCA_000003625.1) using the alignment program hisat2 (14). Non-aligning reads were mapped against the cottontail rabbit papillomavirus reference genome (GenBank Accession NC_001541) by hisat2 to obtain virus-specific reads. To detect viral sequences in human AVM and skin samples, all known pathogenic virus genomes (https://www.ncbi.nlm.nih.gov/genome/viruses/) were combined in a single viral reference genome file and reads were mapped against this reference using Bowtie2 (15). All uniquely mappable reads were re-aligned to the hg19 genome and mapping reads were discarded, yielding only unique reads originating from non-endogenous viral genomes. Read numbers for viral sequences in rabbit and human samples were analyzed quantitatively relative to the total number of reads. For VirCapSeq-VERT analysis, raw data were demultiplexed, Q30-filtered, evaluated by PRINSEQ (v 0.20.2) software (16) and trimmed. Quality-filtered reads were aligned against a host reference database (human genome, including ribosomal RNA and mitochondrial sequences) to remove cellular background and resulting reads de novo assembled using MIRA (v 4.0). Contigs (i.e. assembled overlapping reads yielding a larger segment of the gene) and unique singletons were subjected to homology search using MegaBlast against the GenBank nucleotide database; sequences that showed poor or no homology at the nucleotide level were screened by BLASTX against the viral GenBank protein database. Potential viral sequences from BLASTX analysis were subjected to another round of BLASTX homology search against the entire GenBank protein database to correct for biased e-values and taxonomic misassignments. A positive viral signal was assigned to samples with a read count >10/million quality-filtered, host subtracted reads that distributed to at least three genomic regions.
Results
Prior to processing, representative regions of the samples were FFPE for subsequent use in immunohistochemistry. Figure 1 shows the immunohistochemical staining of the general endothelial cell marker CD31 in different VA to visualize the respective vascular architecture in comparison with normal human skin. RNA-Seq analyses of human AVM and skin control tissues did not support the presence of active viral infection within the tested tissues (Figure 2A). In sharp contrast, viral sequences specific for the cottontail rabbit papilloma virus were readily detected in VX2 tumors of the rabbit, which are known to be transformed and driven by oncoproteins of this papillomavirus, but were absent from normal rabbit skin (Figure 2A). A detailed analysis of the low-level basal reads as seen for skin and AVM (Figure 2A) was performed to evaluate if the putative virus read-levels between skin and AVM were significantly different. For this, read-levels of same viruses instead of viral transcripts in general were compared with each other. Only viruses with more than 10 reads/million (0.001%) were considered for analysis, yielding a total of 11 putative virus candidates. No significant differences were found between putative virus read-levels for skin and AVM (Figure 2B). Furthermore, the presence of human-unrelated putative pathogenic viral reads such as sequences matching pestivirus giraffe-1 (H138) and Saccharomyces cerevisiae killer virus M1 (ScV-M1) underscores these reads as representing unspecific sequences. This is supported by the typical feature of such spurious reads covering only a tiny portion of the respective viral genome, which is illustrated in Figure 2C for the putative reads of murine osteosarcoma virus (MSV).
To account for the presence of low number viral transcripts that might be present in VA tissues and escaped the RNA-Seq approach, we evaluated 10 VA (AVM-2x, HA-1x, LM-1x, VLM-2x, VM-4x) with the highly sensitive virome capture platform VirCapSeq-VERT (10). An average of 36 million raw reads were generated per sample (range=31.4-41 million), yielding a mean of 4 million quality filtered and host subtracted reads (range=3-5.3 million) for contiging (using contigs) and blast analysis.
After bioinformatic assessment, only reads matching viral sequences derived from members of the Retro-, Flavi-, Hepadna- and Hepeviridae families remained at elevated numbers, noted particularly for the family Retroviridae (Figure 3A). Further analysis of the retroviruses revealed sequences for human T-cell lymphotropic virus type 1 (HTLV-1), human immunodeficiency virus type 1 (HIV-1) and human endogenous retrovirus K (HERV-K), respectively, exhibiting substantial numbers particularly for HTLV-1 (Figure 3B). When mapping the putative HTLV-1 reads to the respective genome (17), only one circumscribed region of the genome showed coverage (Figure 3C).
In one sample (VM, HNO-1001-3347, Table I), flaviviral reads did cross the threshold of 0.001% (10 reads/million) and distributed to more than three genomic regions, indicating the presence of a virus (human pegivirus) (Figure 3D). However, the virus is reported to be prevalent in healthy blood donors, not known to cause human pathology and was encountered in only one out of 10 VA specimens (1/4 VM specimens), which does not suggest a role in pathology of vascular anomalies.
Discussion
The pathogenesis of VA is still ill-defined. During the past two decades there were mainly anecdotic reports about gene mutations being associated with development of VA. However, due to the emergence of next-generation sequencing techniques that allow for extensive exome, transcriptome and epigenome analysis, larger groups of patients with VA can now be comprehensively investigated for the presence of mutations that potentially are implicated in development or progression of VA. Among the most prominent candidate genes identified are PIK3CA, rat sarcoma proto oncogene (RAS) and G-protein-coupled receptors (GPCRs) to name just a few (18). Some of the genes or related signaling pathways are druggable and therefore could serve as therapeutic targets.
Apart from the more or less frequently found gene mutations, the initial hit leading to a disturbance of vascular architecture and VA development remains unknown. Many VA manifest shortly after birth or even during pregnancy pointing to a pathogenic factor that is interfering with vascular development during the pregnancy phase. For many years, viral infection was also suspected to play a role in development of VA but for a long time a comprehensive analysis was not possible. This was mainly due to technical limitations allowed investigation of only a limited number of viruses at the same time.
In our study, we used next-generation sequencing-based methods to comprehensively evaluate the presence of viral transcripts in VA tissues. Using RNA-Seq on whole cellular RNA, we did not find a significant number of genuine viral sequences. Similarly, using the highly sensitive VirCapSeq-VERT method that involved the initial capture and enrichment of viral nucleic acids derived from all known vertebrate (including human) pathogenic viruses, we also did not detect any distinct viral transcripts pointing to an association with the tested VA tissues. Borderline levels of pegivirus-specific reads covering more than three genomic regions were detected in one VM but in no other VA. Pegivirus sequences are a common finding in healthy human blood donors and not known to cause human pathology (19). Thus, pegivirus is not considered as being associated with the pathogenesis of VA. Our data therefore does not support the presence of an active viral infection as being a driver of already developed VA. However, since the investigations (except for samples with ID HNO-1001-1756 and -2315) considered VA tissues from adults (average age 38 years; standard deviation ±22 years; n=13; see Table I), it cannot be excluded that at an earlier stage, e.g. during pregnancy, a viral infection had taken place that resulted in a disturbance of vascular development and possibly even in the gene mutations found associated with VA (18). There are several examples in the literature in which a viral infection during pregnancy affected the vascular system. Infection with herpesvirus type 6, for example, affects vessels of the fetoplacental system (20).
Types of VA tissues used for viral analysis. Representative microscopical images are shown of immunohistological staining using an anti CD31 specific antibody to visualize CD31+ (brown precipitate) endothelial cells. AVM: Arteriovenous malformation, HA: hemangioma, LM: lymphatic malformation, VM: venous malformation, VLM: venous-lymphatic malformation. Tissues used for CD31 immunohistochemistry are HNO-1000-720, -2590, -23015, -601, -2873, -3016 (see Table I).
RNA-Seq analysis does not support the presence of viral transcripts in tissues from arteriovenous malformations (AVM). A: RNA-Seq analysis in human healthy skin, AVM tissues, and rabbit (rb) VX2 tumor tissues, and skin. No significant viral transcript levels were detected in human tissues. Cottontail rabbit papillomavirus (CRPV)-specific transcripts were readily detectable in VX2 tumor tissues (positive control) but not in rabbit skin (negative control). B: Detailed evaluation of the seeming basal viral reads as seen for skin and AVM tissues in (A) did not show statistical significant differences between the two groups. C: Genome delineation adapted from the graphical view created by the nucleotide database of the National Center for Biotechnology Information (NCBI, Rockville Pike, MD, USA). Only limited read coverage is seen for the murine osteosarcoma virus (MSV) genome, pointing to the presence of unspecific virus unrelated reads. BBMV (RNA 3): Broad bean mottle virus RNA 3; FSV: Fujinami sarcoma virus); ScV-M1: Saccharomyces cerevisiae killer virus M1; HLSV: Hibiscus latent Singapore virus; WMSV: wooly monkey sarcoma virus; EMCV: encephalomyocarditis virus; FHMPV: fathead minnow picornavirus; A-MuLV: Abelson murine leukemia virus; MSV: murine osteosarcoma virus.
VirCapSeq-VERT. Analysis did not reveal a distinct association of viral transcripts with VA. After bioinformatic filtration (A), only a limited number of viral reads, mostly for viruses of the family Retroviridae, remained detectable. Of the retroviruses involved, human T-cell lymphotropic virus type 1 (HTLV-1) exhibited a substantial number of reads (B), which however covered only a tiny portion of the proviral genome (C). Threshold level pegivirus-specific reads (D) covered more than three genomic regions in one VM. AVM: Arteriovenous malformation, HA: hemangioma, LM: lymphatic malformation, VM: venous malformation, VLM: venous-lymphatic malformation; HIV-1: human immunodeficiency virus type 1; HERV-K: human endogenous retrovirus K.
Similarly, cytomegalovirus (CMV) infection results in a disturbance of placental angiogenesis (21), which was related to infection of placental pericytes (22). Furthermore, endothelial cells infected with CMV displayed an aberrant vascular morphology (23) and exhibited elevated levels of leukotriene B4 and 5-lipoxygenase (24). In addition, CMV infection was found to inhibit lymphangiogenesis and vascular remodeling in a human placenta model (25) and resulted in impaired vascular function in a mouse model (26). Interestingly, varicella zoster virus infection was suspected to play a role in the development of tufted angiomas at a higher age (27).
In conclusion, our data do not support the presence of an active viral infection in established VA of the head and neck region. However, a viral infection occurring during pregnancy might still be the initial pathogenic trigger that induces a disturbance of endothelial–vascular development leading to the appearance of VA. Evaluation of this aspect will be a task for future studies.
Acknowledgements
Technical assistance by Ms. R. Peldszus and Ms. M. Sadowski (Department of Otorhinolaryngology, Head and Neck Surgery, University Hospital Giessen and Marburg, Campus Marburg, Marburg, Germany) was greatly appreciated. The Authors thank the team of Dr. Vladimir Beneš (Genomics Core Facility, European Molecular Biology Laboratory, Heidelberg, Germany) for their support during RNA-Seq. R.M. highly appreciates the discussions he had with K. Jain, N. Bhuva, C. Guo, A. Caciula (Columbia University, New York City, NY, USA) during his visit of the Lipkin Laboratory. Parts of the study were supported by the Kempkes Stiftung (Philipps-Universität Marburg, Marburg, Germany)
Footnotes
- Received August 31, 2018.
- Revision received September 20, 2018.
- Accepted September 28, 2018.
- Copyright© 2018, International Institute of Anticancer Research (Dr. George J. Delinasios), All rights reserved