Contents

1 Installation

EventPointer can be installed from Bioconductor using the following lines:

source("http://www.bioconductor.org/biocLite.R")
biocLite("EventPointer")

2 Introduction

EventPointer R package provides users a simplified way to identify, classify and visualize alternative splicing events. The steps required by the algorithm are almost identical for both technologies. The algorithm only differs in its inital step.

Figure 1. Definition of alternative splicing events for EventPointer

Figure 1. Definition of alternative splicing events for EventPointer

Figure 2 shows each of the main steps in a graphical layout.

This vignette is divided in two sections. In the first one, the complete analysis for junction arrays is described and the second one describes the analysis for RNA-Seq data.

To cite EventPointer:

Figure 2. Overview of EventPointer

Figure 2. Overview of EventPointer

3 Analysis of junction arrays

3.1 Overview of junction arrays

EventPointer is prepared to work with different Affymetrix arrays, such as: HTA 2.0, Clariom-D, RTA and MTA.

To build the CDF file (to work under the aroma.affymetrix framework), EventPointer requires a GTF file with the reference transcriptome information. In case not provided, the algorithm automatically downloads the required information from different databases such as ENSEMBL or UCSC. As the probes for the HTA 2.0 array are mapped to the HG19 genome, the latests versions of the ensembl and ucsc genome, mapped to the hg19 version, will be downloaded. For the other arrays, the following genomes are used: ClariomD = GRCh38, MTA = mm10 and RTA = rn6.

The required files are:

  1. Exon probes genomic position (Tab separated .txt file)
  2. Junction probes genomic position (Tab separated .txt file)
  3. Reference transcriptome (GTF file)

Files 1 and 2 contain probe information for the array. These files and the corresponding CDF files can be downloaded from the following URL: EventPointer Dropbox

The format of these files is briefly explained in the following paragraphs:

For the Exon Probes, the file corresponds to a tab separated .txt file composed of 11 columns that include: Probe Id, X coordinate in the array, Y coordinate in the array, Transcript Cluster Id, Probeset Id, Probeset name, Probe sequence, chromosome, start, end and strand.

The Junction probes file is also a tab separated .txt composed of 10 columns: Probe Id, X coordinate in the array, Y coordinate in the array, Transcript Cluster Id, Probeset Id, Probeset name, Probe sequence, chromosome and probe alignments.

For the GTF the standard format is used. (For more information see GTF)

This vignette uses the probes and annotation file for the DONSON gene in order to run the examples in a short amount of time. To generate the corresponding CDF file for the whole genome, the files from the Dropbox link must be used.

Note: It is advisable to work with reference GTF files, as the probes are annotated to them. However, if other database is used, the algorithm will only include probes that are mapped to such database.

3.2 CDF file creation

This step can be skipped if the corresponding CDF file is doownloaded from the Dropbox link. The available CDF files were generated using the GTF files for each of the arrays, if users want to generate a CDF for other databases (Ensembl or UCSC), this step should be used.

The function CDFfromGTF generates the CDF file used afterwards in the aroma.affymetrix pre-processing pipeline. Internally, it calls flat2cdf function written by Elizabeth Purdom. More information about it can be seen in the following webpage: Create CDF from scratch

The required input for the function is described below:

  • input : Reference transcriptome. Available options are : “Ensembl”,“UCSC” , “AffyGTF” and “CustomGTF”.
  • inputFile: If input is “AffyGTF” or “CustomGTF”, inputFile should point to the GTF file to be used.
  • PSR: Path to the Exon probes txt file
  • Junc: Path to the Junction probes txt file
  • PathCDF: Directory where the output will be saved
  • microarray: Microarray used to create the CDF file. Must be one of: “HTA-2_0”, “ClariomD”, “RTA” or “MTA”.

This function takes a couple of hours to complete (depending on the computer), and creates the following files:

  1. EventsFound.txt : Tab separated file with all the information of all the alternative splcing events found.
  2. .flat file : Used to build the corresponding CDF file.
  3. .CDF file: Output required for the aroma.affymetrix preprocessing pipeline.

The following code chunks show examples on how to run the function using part of the Affymetrix GTF and the example data included in the package. This data corresponds to the Affymetrix HTA 2.0 GTF representing only the DONSON gene and the probes that are mapped to it.

Using Affymetrix GTF as reference transcriptome


# Set input variables
   PathFiles<-system.file("extdata",package="EventPointer")
   DONSON_GTF<-paste(PathFiles,"/DONSON.gtf",sep="")
   PSRProbes<-paste(PathFiles,"/PSR_Probes.txt",sep="")
   JunctionProbes<-paste(PathFiles,"/Junction_Probes.txt",sep="")
   Directory<-tempdir()
   array<-"HTA-2_0"

# Run the function

   CDFfromGTF(input="AffyGTF",inputFile=DONSON_GTF,
              PSR=PSRProbes,Junc=JunctionProbes,
              PathCDF=Directory,microarray=array)
## Creating SG Information...
## Import genomic features from the file as a GRanges object ... OK
## Prepare the 'metadata' data frame ... OK
## Make the TxDb object ... OK
## 
## Reading Information On Probes...Done
## Indexing Genes and Probes...Done
## Finding Events...
## 
  |                                                                         
  |                                                                   |   0%
  |                                                                         
  |===================================================================| 100%
## 
## Creating .flat ...DoneReading TXT file ... Done.
## Splitting TXT file into units ... split into 1 initial chunks ... unwrapped into 12 chunks ... Done.
## Creating structure for 12 units (dot=250):
## 
## Writing CDF file...
##   Pathname: HTA-2_0.cdf
##   Writes CDF header and unit names...
##   Writes CDF header and unit names...done
##   Writes QC units...
##     Units left: 0
##   Writes QC units...done
##   Writes 12 units...
##     Units left: 0
##   Writes 12 units...done
## Writing CDF file...done

Note: Both the .flat and .CDF file take large amounts of space in the hard drive, it is recommended that the directory has at least 1.5 GB of free space.

Figure 3 shows a screen shot with the output of the CDFfromGTF function

Figure 3. Output of CDFfromGTF

Figure 3. Output of CDFfromGTF

Once the file is created, the name of the cdf file can be changed. Internally, the algorithm gives the name as HTA-2_0, but the actual name of the file can be different. In the rest of the vignette, we have renamed our CDF file as EP_HTA-2_0.

Once the CDF file has been created, it can be used for the analysis of different experiments.

3.3 Statistical Analysis

3.3.1 aroma.affymetrix pre-processing pipeline

For microarray data, a pre-processing pipeline must be applied to the .cel files of the experiment. Usually this pre-processing is done using the aroma.affymetrix R package. This procedure normalizes and summarizes the expression of the different probesets into single values.

The aroma.affymetrix R package provides users different functions to work with affymetrix arrays. The functions are used to extract all the information contained in the .cel files and to do all the required pre-processing steps such as background correction, normalization and summarization. The package requires a certain folder structure in order to work correctly. For more information about aroma.affymetrix visit the webpage:aroma project

The following code chunk displays the pipeline used to get the results required by EventPointer after the pre-processing using aroma.affymetrix. The code is not intended to be a runable example, but just to show users the settings and functions that can be used. In order for users to have a runable example, the corrrespoding folder structure and files are required.

# Simple example of Aroma.Affymetrix Preprocessing Pipeline

verbose <- Arguments$getVerbose(-8);
timestampOn(verbose);
projectName <- "Experiment"
cdfGFile <- "EP_HTA-2_0,r"
cdfG <- AffymetrixCdfFile$byChipType(cdfGFile)
cs <- AffymetrixCelSet$byName(projectName, cdf=cdfG)
bc <- NormExpBackgroundCorrection(cs, method="mle", tag=c("*","r11"));
csBC <- process(bc,verbose=verbose,ram=20);
qn <- QuantileNormalization(csBC, typesToUpdate="pm");
csN <- process(qn,verbose=verbose,ram=20);
plmEx <- ExonRmaPlm(csN, mergeGroups=FALSE)
fit(plmEx, verbose=verbose)
cesEx <- getChipEffectSet(plmEx)
ExFit <- extractDataFrame(cesEx, addNames = TRUE)

3.3.2 EventPointer function

EventPointer is the main function of the algorithm

The function requires the following parameters:

  • Design : The design matrix for the experiment.
  • Contrast: The contrast matrix for the experiment.
  • ExFit: [aroma.affymetrix] pre-processed variable after using extractDataFrame(affy, addNames=TRUE)
  • Eventstxt: Path to the EventsFound.txt file generated by CDFfromGTF function.
  • Filter: Boolean variable to indicate if an expression filter is applied.
  • Qn: Quantile used to filter the events (Bounded between 0-1, Q1 would be 0.25).
  • Statistic: Statistical test to identify differential splicing events, must be one of : “LogFC”,“Dif_LogFC” and “DRS”.
  • PSI: Boolean variable to indicate if \(\Delta \Psi\) should be calculated for every splicing event.

If the Filter variable is TRUE, the following is performed:

The summarized expression value of all the reference paths is obtained and the threshold is set depending on the Qn value used.

An event is considered if at least one sample in all paths is expressed above the threshold.

The rest of the events are not shown unless the Filter variable is set to FALSE

   data(ArraysData)

   Dmatrix<-matrix(c(1,1,1,1,0,0,1,1),nrow=4,ncol=2,byrow=FALSE)
   Cmatrix<-t(t(c(0,1)))
   EventsFound<-paste(system.file("extdata",package="EventPointer"),"/EventsFound.txt",sep="")
   
   Events<-EventPointer(Design=Dmatrix,
                      Contrast=Cmatrix,
                      ExFit=ArraysData,
                      Eventstxt=EventsFound,
                      Filter=FALSE,
                      Qn=0.25,
                      Statistic="LogFC",
                      PSI=TRUE)
## 09:17:24 PM Running EventPointer:  
##  ** Statistical Analysis: Logarithm of the fold change of both isoforms 
##  ** Delta PSI will be calculated 
##  ** No expression filter 
##             ---------------------------------------------------------------- 
##  ** Calculating PSI...Done 
##  ** Running Statistical Analysis...Done 
## 
##  09:17:25 PM  Analysis Completed

Table 1 displays the output of EventPointer function


(#tab:EP_Arrays_Res_Table)Table 1: EventPointer Arrays results
Gene name Event Type Genomic Position Splicing Z Value Splicing Pvalue Delta PSI
TC21001058.hg_9 TC21001058.hg Complex Event 21:34951868-34956896 7.16865 0 -0.58200
TC21001058.hg_2 TC21001058.hg Complex Event 21:34954361-34956896 6.82303 0 0.00000
TC21001058.hg_12 TC21001058.hg Complex Event 21:34958487-34960627 6.33509 0 -0.09482
TC21001058.hg_3 TC21001058.hg Complex Event 21:34954552-34956896 6.11308 0 -0.45633
TC21001058.hg_6 TC21001058.hg Complex Event 21:34958487-34960627 -6.08339 0 0.43586

The columns of the data.frame are:

  • Gene name : Gene identifier
  • Event Type: Type of alternative splicing event (Cassette, Alternative 3’,Alternative 5’, Alternative Last Exon, Alternative First Exon, Mutually Exclusive Exons or Complex Event)
  • Genomic Position: Genomic Coordinates for the event
  • Splicing Z Value: Corresponding Z value for the statistical test performed
  • Splicing P Value: Corresponding P-value for the statistical test performed
  • Delta PSI: \(\Delta \Psi\) parameter for each event

3.4 IGV visualization

EventPointer creates two different GTF files to visualize the alternative splicing events. Figure 5 displays the cassette exon for the COPS7A gene (4th ranked event in Table 1). In the IGV visualization, the reference path is shown in gray color, the path 1 in red and path 2 in green. Below the transcripts, the different probes that are measuring each of the paths are represented in the same colors.

Figure 5. EventPointer visualization using IGV

Figure 5. EventPointer visualization using IGV

To create the GTF files, the algorithm uses the EventPointer_IGV functions with the following parameters:

  • Events : Data.frame generated by EventPointer with the events to be included in the GTF file.
  • input: Reference transcriprome. Must be one of: “Ensembl”, “UCSC” , “AffyGTF” or “CustomGTF”.
  • inputFile: If input is “AffyGTF” or “CustomGTF”, inputFile should point to the GTF file to be used.
  • PSR: Path to the Exon probes txt file.
  • Junc: Path to the Junction probes txt file.
  • PathGTF: Directory where to write the GTF files.
  • EventsFile: Path to EventsFound.txt file generated with CDFfromGTF function.
  • microarray: Microarray used to create the CDF file. Must be one of: HTA-2_0, ClariomD, RTA or MTA

The inital process of the function is slow, as the s