Contents

1 Background

The data used in this workflow comes from an RNA-seq experiment where airway smooth muscle cells were treated with dexamethasone, a synthetic glucocorticoid steroid with anti-inflammatory effects (Himes et al. 2014). Glucocorticoids are used, for example, by people with asthma to reduce inflammation of the airways. In the experiment, four human airway smooth muscle cell lines were treated with 1 micromolar dexamethasone for 18 hours. For each of the four cell lines, we have a treated and an untreated sample. For more description of the experiment see the article, PubMed entry 24926665, and for raw data see the GEO entry GSE52778.

2 Data

FastQ files with a small subset of the reads can be found on https://github.com/statOmics/SGA2019/tree/data-rnaseq

library("Rsubread")
## Warning: package 'Rsubread' was built under R version 3.6.1
library("GEOquery")
## Loading required package: Biobase
## Loading required package: BiocGenerics
## Loading required package: parallel
## 
## Attaching package: 'BiocGenerics'
## The following objects are masked from 'package:parallel':
## 
##     clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
##     clusterExport, clusterMap, parApply, parCapply, parLapply,
##     parLapplyLB, parRapply, parSapply, parSapplyLB
## The following objects are masked from 'package:stats':
## 
##     IQR, mad, sd, var, xtabs
## The following objects are masked from 'package:base':
## 
##     Filter, Find, Map, Position, Reduce, anyDuplicated, append,
##     as.data.frame, basename, cbind, colnames, dirname, do.call,
##     duplicated, eval, evalq, get, grep, grepl, intersect,
##     is.unsorted, lapply, mapply, match, mget, order, paste, pmax,
##     pmax.int, pmin, pmin.int, rank, rbind, rownames, sapply,
##     setdiff, sort, table, tapply, union, unique, unsplit, which,
##     which.max, which.min
## Welcome to Bioconductor
## 
##     Vignettes contain introductory material; view with
##     'browseVignettes()'. To cite Bioconductor, see
##     'citation("Biobase")', and for packages 'citation("pkgname")'.
## Setting options('download.file.method.GEOquery'='auto')
## Setting options('GEOquery.inmemory.gpl'=FALSE)
library("tidyverse")
## -- Attaching packages -------------------- tidyverse 1.2.1 --
## v ggplot2 3.2.1     v purrr   0.3.2
## v tibble  2.1.3     v dplyr   0.8.2
## v tidyr   0.8.3     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.4.0
## -- Conflicts ----------------------- tidyverse_conflicts() --
## x ggplot2::Position() masks BiocGenerics::Position(), base::Position()
## x dplyr::combine()    masks Biobase::combine(), BiocGenerics::combine()
## x dplyr::filter()     masks stats::filter()
## x dplyr::lag()        masks stats::lag()
getwd()
## [1] "/Users/lclement/Downloads/airwayLecture"

We will use the Rsubread read mapper because that is avaible in R for all platforms (Linux, Windows and Mac). For real projects I prefer the use of STAR.

3 Get info on experiment

3.1 Get info on samples

Get all info from GEO. get sample info via getGEO (info from samples)

gse<-getGEO("GSE52778")
## Found 1 file(s)
## GSE52778_series_matrix.txt.gz
## Parsed with column specification:
## cols(
##   ID_REF = col_character(),
##   GSM1275862 = col_character(),
##   GSM1275863 = col_character(),
##   GSM1275864 = col_character(),
##   GSM1275865 = col_character(),
##   GSM1275866 = col_character(),
##   GSM1275867 = col_character(),
##   GSM1275868 = col_character(),
##   GSM1275869 = col_character(),
##   GSM1275870 = col_character(),
##   GSM1275871 = col_character(),
##   GSM1275872 = col_character(),
##   GSM1275873 = col_character(),
##   GSM1275874 = col_character(),
##   GSM1275875 = col_character(),
##   GSM1275876 = col_character(),
##   GSM1275877 = col_character()
## )
## File stored at:
## /var/folders/p1/3js9hvbs473g1klcmqm0d8wm0000gn/T//RtmpA9fc49/GPL11154.soft
length(gse)
## [1] 1

There is one object for this experiment.

pdata<-pData(gse[[1]])
pdata$SampleName<-rownames(pdata) %>% as.factor
head(pdata)[1:6,]
##                        title geo_accession                status
## GSM1275862  N61311_untreated    GSM1275862 Public on Jan 01 2014
## GSM1275863        N61311_Dex    GSM1275863 Public on Jan 01 2014
## GSM1275864        N61311_Alb    GSM1275864 Public on Jan 01 2014
## GSM1275865    N61311_Alb_Dex    GSM1275865 Public on Jan 01 2014
## GSM1275866 N052611_untreated    GSM1275866 Public on Jan 01 2014
## GSM1275867       N052611_Dex    GSM1275867 Public on Jan 01 2014
##            submission_date last_update_date type channel_count
## GSM1275862     Nov 26 2013      May 15 2019  SRA             1
## GSM1275863     Nov 26 2013      May 15 2019  SRA             1
## GSM1275864     Nov 26 2013      May 15 2019  SRA             1
## GSM1275865     Nov 26 2013      May 15 2019  SRA             1
## GSM1275866     Nov 26 2013      May 15 2019  SRA             1
## GSM1275867     Nov 26 2013      May 15 2019  SRA             1
##                       source_name_ch1 organism_ch1
## GSM1275862 airway smooth muscle cells Homo sapiens
## GSM1275863 airway smooth muscle cells Homo sapiens
## GSM1275864 airway smooth muscle cells Homo sapiens
## GSM1275865 airway smooth muscle cells Homo sapiens
## GSM1275866 airway smooth muscle cells Homo sapiens
## GSM1275867 airway smooth muscle cells Homo sapiens
##                           characteristics_ch1
## GSM1275862               treatment: Untreated
## GSM1275863           treatment: Dexamethasone
## GSM1275864               treatment: Albuterol
## GSM1275865 treatment: Albuterol_Dexamethasone
## GSM1275866               treatment: Untreated
## GSM1275867           treatment: Dexamethasone
##                               characteristics_ch1.1 characteristics_ch1.2
## GSM1275862 tissue: human airway smooth muscle cells           ercc_mix: -
## GSM1275863 tissue: human airway smooth muscle cells           ercc_mix: -
## GSM1275864 tissue: human airway smooth muscle cells           ercc_mix: -
## GSM1275865 tissue: human airway smooth muscle cells           ercc_mix: -
## GSM1275866 tissue: human airway smooth muscle cells           ercc_mix: -
## GSM1275867 tissue: human airway smooth muscle cells           ercc_mix: 1
##            characteristics_ch1.3                 characteristics_ch1.4
## GSM1275862     cell line: N61311 cell type: airway smooth muscle cells
## GSM1275863     cell line: N61311 cell type: airway smooth muscle cells
## GSM1275864     cell line: N61311 cell type: airway smooth muscle cells
## GSM1275865     cell line: N61311 cell type: airway smooth muscle cells
## GSM1275866    cell line: N052611 cell type: airway smooth muscle cells
## GSM1275867    cell line: N052611 cell type: airway smooth muscle cells
##                                                                                                                                               treatment_protocol_ch1
## GSM1275862 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
## GSM1275863 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
## GSM1275864 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
## GSM1275865 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
## GSM1275866 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
## GSM1275867 Cells from each donor were treated with 1 \316\274M dexamethasone,  1 \316\274M albuterol, both dexamethasone and albuterol, or control vehicle for 18 h.
##                                                                                                                                                                                                                        growth_protocol_ch1
## GSM1275862 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
## GSM1275863 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
## GSM1275864 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
## GSM1275865 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
## GSM1275866 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
## GSM1275867 Primary ASM cells were isolated from four white aborted lung transplant donors with no chronic illness. Passages 4\302\240to 7 ASM cells maintained in Ham's F12 medium supplemented with 10% FBS were used in all experiments.
##            molecule_ch1
## GSM1275862    total RNA
## GSM1275863    total RNA
## GSM1275864    total RNA
## GSM1275865    total RNA
## GSM1275866    total RNA
## GSM1275867    total RNA
##                                                                                     extract_protocol_ch1
## GSM1275862 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
## GSM1275863 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
## GSM1275864 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
## GSM1275865 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
## GSM1275866 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
## GSM1275867 Total RNA was extracted using the miRNAeasy mini kit (Qiagen Sciences, Inc., Germantown, MD).
##                                                    extract_protocol_ch1.1
## GSM1275862 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
## GSM1275863 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
## GSM1275864 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
## GSM1275865 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
## GSM1275866 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
## GSM1275867 TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA).
##            taxid_ch1                                   data_processing
## GSM1275862      9606 Illumina Casava1.8 software used for basecalling.
## GSM1275863      9606 Illumina Casava1.8 software used for basecalling.
## GSM1275864      9606 Illumina Casava1.8 software used for basecalling.
## GSM1275865      9606 Illumina Casava1.8 software used for basecalling.
## GSM1275866      9606 Illumina Casava1.8 software used for basecalling.
## GSM1275867      9606 Illumina Casava1.8 software used for basecalling.
##                                                                                                                                                              data_processing.1
## GSM1275862 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
## GSM1275863 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
## GSM1275864 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
## GSM1275865 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
## GSM1275866 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
## GSM1275867 Sequenced reads were trimmed by 12bp then mapped to hg19 whole genome using tophat v.2.0.4 with parameters -G [hg19_ERCC_gtf] --no-novel-juncs --transcriptome-only
##                                                                                                                                                                                                    data_processing.2
## GSM1275862 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
## GSM1275863 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
## GSM1275864 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
## GSM1275865 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
## GSM1275866 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
## GSM1275867 For each sample, Cufflinks (v.2.0.2) was used to quantify External RNA Controls Consortium (ERCC) Spike-In and hg19 transcripts based on reads that mapped to the provided hg19 and ERCC reference files.
##                                                                                                                                                                                                                                                                                                                                                                                                                   data_processing.3
## GSM1275862 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
## GSM1275863 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
## GSM1275864 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
## GSM1275865 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
## GSM1275866 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
## GSM1275867 Differential expression of genes and transcripts in samples treated with Dex vs. Untreated samples was obtained using Cuffdiff (v.2.0.2) with the quantified transcripts computed by Cufflinks (v.2.0.2), while applying bias correction. Options: -p 12 --min-outlier-p 0.05 -b [hg19_reference_index] -u [hg19_gtf] -L Dex,Control  [List of accepted_hits.bam for each sample sorted by condition (i.e. Dex,Control)]
##                                                                       data_processing.4
## GSM1275862 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
## GSM1275863 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
## GSM1275864 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
## GSM1275865 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
## GSM1275866 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
## GSM1275867 CummeRbund R package (v.0.1.3) was used to create FPKM matrix of all samples
##             data_processing.5
## GSM1275862 Genome_build: hg19
## GSM1275863 Genome_build: hg19
## GSM1275864 Genome_build: hg19
## GSM1275865 Genome_build: hg19
## GSM1275866 Genome_build: hg19
## GSM1275867 Genome_build: hg19
##                                                                                                                                                                                                                                                             data_processing.6
## GSM1275862 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
## GSM1275863 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
## GSM1275864 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
## GSM1275865 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
## GSM1275866 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
## GSM1275867 Supplementary_files_format_and_content: tab-delimited text files on Series record include FPKM values for each sample (ALL_Sample_FPKM_Matrix.txt) and gene-based differential expression results for Dex vs. Untreated condition (Dex_vs_Untreated_gene_exp.diff)
##            platform_id   contact_name                 contact_email
## GSM1275862    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
## GSM1275863    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
## GSM1275864    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
## GSM1275865    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
## GSM1275866    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
## GSM1275867    GPL11154 Blanca,E,Himes bhimes@pennmedicine.upenn.edu
##            contact_laboratory
## GSM1275862  402 Blockley Hall
## GSM1275863  402 Blockley Hall
## GSM1275864  402 Blockley Hall
## GSM1275865  402 Blockley Hall
## GSM1275866  402 Blockley Hall
## GSM1275867  402 Blockley Hall
##                                                   contact_department
## GSM1275862 Department of Biostatistics, Epidemiology and Informatics
## GSM1275863 Department of Biostatistics, Epidemiology and Informatics
## GSM1275864 Department of Biostatistics, Epidemiology and Informatics
## GSM1275865 Department of Biostatistics, Epidemiology and Informatics
## GSM1275866 Department of Biostatistics, Epidemiology and Informatics
## GSM1275867 Department of Biostatistics, Epidemiology and Informatics
##                     contact_institute contact_address contact_city
## GSM1275862 University of Pennsylvania 423 Guardian Dr Philadelphia
## GSM1275863 University of Pennsylvania 423 Guardian Dr Philadelphia
## GSM1275864 University of Pennsylvania 423 Guardian Dr Philadelphia
## GSM1275865 University of Pennsylvania 423 Guardian Dr Philadelphia
## GSM1275866 University of Pennsylvania 423 Guardian Dr Philadelphia
## GSM1275867 University of Pennsylvania 423 Guardian Dr Philadelphia
##            contact_state contact_zip/postal_code contact_country
## GSM1275862            PA                   19104             USA
## GSM1275863            PA                   19104             USA
## GSM1275864            PA                   19104             USA
## GSM1275865            PA                   19104             USA
## GSM1275866            PA                   19104             USA
## GSM1275867            PA                   19104             USA
##            data_row_count    instrument_model library_selection
## GSM1275862              0 Illumina HiSeq 2000              cDNA
## GSM1275863              0 Illumina HiSeq 2000              cDNA
## GSM1275864              0 Illumina HiSeq 2000              cDNA
## GSM1275865              0 Illumina HiSeq 2000              cDNA
## GSM1275866              0 Illumina HiSeq 2000              cDNA
## GSM1275867              0 Illumina HiSeq 2000              cDNA
##            library_source library_strategy
## GSM1275862 transcriptomic          RNA-Seq
## GSM1275863 transcriptomic          RNA-Seq
## GSM1275864 transcriptomic          RNA-Seq
## GSM1275865 transcriptomic          RNA-Seq
## GSM1275866 transcriptomic          RNA-Seq
## GSM1275867 transcriptomic          RNA-Seq
##                                                                  relation
## GSM1275862 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422669
## GSM1275863 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422675
## GSM1275864 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422668
## GSM1275865 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422667
## GSM1275866 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422678
## GSM1275867 BioSample: https://www.ncbi.nlm.nih.gov/biosample/SAMN02422670
##                                                      relation.1
## GSM1275862 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384345
## GSM1275863 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384346
## GSM1275864 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384347
## GSM1275865 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384348
## GSM1275866 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384349
## GSM1275867 SRA: https://www.ncbi.nlm.nih.gov/sra?term=SRX384350
##            supplementary_file_1 cell line:ch1              cell type:ch1
## GSM1275862                 NONE        N61311 airway smooth muscle cells
## GSM1275863                 NONE        N61311 airway smooth muscle cells
## GSM1275864                 NONE        N61311 airway smooth muscle cells
## GSM1275865                 NONE        N61311 airway smooth muscle cells
## GSM1275866                 NONE       N052611 airway smooth muscle cells
## GSM1275867                 NONE       N052611 airway smooth muscle cells
##            ercc_mix:ch1                       tissue:ch1
## GSM1275862            - human airway smooth muscle cells
## GSM1275863            - human airway smooth muscle cells
## GSM1275864            - human airway smooth muscle cells
## GSM1275865            - human airway smooth muscle cells
## GSM1275866            - human airway smooth muscle cells
## GSM1275867            1 human airway smooth muscle cells
##                      treatment:ch1 SampleName
## GSM1275862               Untreated GSM1275862
## GSM1275863           Dexamethasone GSM1275863
## GSM1275864               Albuterol GSM1275864
## GSM1275865 Albuterol_Dexamethasone GSM1275865
## GSM1275866               Untreated GSM1275866
## GSM1275867           Dexamethasone GSM1275867

3.2 Get info on sequencing files

Download SRA info. To link sample info to info sequencing: Go to corresponding SRA page and save the information via the “Send to: File button” This file can also be used to make a script to download sequencing files from the web. Note that sra files can be converted to fastq files via the fastq-dump function of the sra-tools.

sraInfo<-read.csv("SraRunInfo.csv")
levels(pdata$SampleName)==levels(sraInfo$SampleName)
##  [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [15] TRUE TRUE

SampleNames are can be linked.

pdata<-merge(pdata,sraInfo,by="SampleName")
pdata$Run
##  [1] SRR1039508 SRR1039509 SRR1039510 SRR1039511 SRR1039512 SRR1039513
##  [7] SRR1039514 SRR1039515 SRR1039516 SRR1039517 SRR1039518 SRR1039519
## [13] SRR1039520 SRR1039521 SRR1039522 SRR1039523
## 16 Levels: SRR1039508 SRR1039509 SRR1039510 SRR1039511 ... SRR1039523

The run is also the name of the SRA file so we will be able to link alignment file name to the experiment via the SRA file info.

4 Path to Genome index

All reads in the subsampled fastq files map to chromosome 1. We have created an index for the chromosome 1 of the Human reference genome in an other R-markdown file. We will reuse it here.

indexName<-"airway_index/homo_sapiens_GRCh38_dna_chromosome_1_rsubread"

4.1 Set path to reads and output

fls1 <- list.files("fastq", full=TRUE,pattern="1.fastq")
fls2 <- list.files("fastq", full=TRUE,pattern="2.fastq")
system("mkdir bamDir")
bamDir<-"./bamDir"
bamfls<-paste0(bamDir,"/",substr(basename(fls1),1,10),".bam")

4.2 Readmapping

The offset for the phred scores is 64. We find info on illumina incoding in quality control step of fastQC.

phredOffset<-33
align(index=indexName,readfile1=fls1,readfile2=fls2,input_format="FASTQ",output_format="BAM",output_file=bamfls,phredOffset=phredOffset)
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039508_subset_1.fastq                                  ||
## || Input file 2  : SRR1039508_subset_2.fastq                                  ||
## || Output file   : SRR1039508.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:15:52, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 139 bp                                       ||
## ||   90% completed, 0.1 mins elapsed, rate=1.0k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 7097                                         ||
## ||                      Mapped : 7090 (99.9%)                                 ||
## ||             Uniquely mapped : 7086                                         ||
## ||               Multi-mapping : 4                                            ||
## ||                                                                            ||
## ||                    Unmapped : 7                                            ||
## ||                                                                            ||
## ||             Properly paired : 5248                                         ||
## ||         Not properly paired : 1842                                         ||
## ||                   Singleton : 45                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 9                                            ||
## ||  Unexpected fragment length : 1784                                         ||
## ||       Unexpected read order : 4                                            ||
## ||                                                                            ||
## ||                      Indels : 38                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039509_subset_1.fastq                                  ||
## || Input file 2  : SRR1039509_subset_2.fastq                                  ||
## || Output file   : SRR1039509.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:15:59, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 149 bp                                       ||
## ||   89% completed, 0.1 mins elapsed, rate=1.7k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 7168                                         ||
## ||                      Mapped : 7157 (99.8%)                                 ||
## ||             Uniquely mapped : 7153                                         ||
## ||               Multi-mapping : 4                                            ||
## ||                                                                            ||
## ||                    Unmapped : 11                                           ||
## ||                                                                            ||
## ||             Properly paired : 5123                                         ||
## ||         Not properly paired : 2034                                         ||
## ||                   Singleton : 55                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 15                                           ||
## ||  Unexpected fragment length : 1961                                         ||
## ||       Unexpected read order : 3                                            ||
## ||                                                                            ||
## ||                      Indels : 28                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039512_subset_1.fastq                                  ||
## || Input file 2  : SRR1039512_subset_2.fastq                                  ||
## || Output file   : SRR1039512.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:03, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 136 bp                                       ||
## ||   86% completed, 0.1 mins elapsed, rate=2.1k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 8478                                         ||
## ||                      Mapped : 8461 (99.8%)                                 ||
## ||             Uniquely mapped : 8457                                         ||
## ||               Multi-mapping : 4                                            ||
## ||                                                                            ||
## ||                    Unmapped : 17                                           ||
## ||                                                                            ||
## ||             Properly paired : 6187                                         ||
## ||         Not properly paired : 2274                                         ||
## ||                   Singleton : 58                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 8                                            ||
## ||  Unexpected fragment length : 2202                                         ||
## ||       Unexpected read order : 6                                            ||
## ||                                                                            ||
## ||                      Indels : 34                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039513_subset_1.fastq                                  ||
## || Input file 2  : SRR1039513_subset_2.fastq                                  ||
## || Output file   : SRR1039513.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:06, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 147 bp                                       ||
## ||   88% completed, 0.1 mins elapsed, rate=1.9k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 7514                                         ||
## ||                      Mapped : 7507 (99.9%)                                 ||
## ||             Uniquely mapped : 7503                                         ||
## ||               Multi-mapping : 4                                            ||
## ||                                                                            ||
## ||                    Unmapped : 7                                            ||
## ||                                                                            ||
## ||             Properly paired : 5407                                         ||
## ||         Not properly paired : 2100                                         ||
## ||                   Singleton : 76                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 10                                           ||
## ||  Unexpected fragment length : 2005                                         ||
## ||       Unexpected read order : 9                                            ||
## ||                                                                            ||
## ||                      Indels : 25                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039516_subset_1.fastq                                  ||
## || Input file 2  : SRR1039516_subset_2.fastq                                  ||
## || Output file   : SRR1039516.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:10, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 139 bp                                       ||
## ||   85% completed, 0.1 mins elapsed, rate=2.1k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 8771                                         ||
## ||                      Mapped : 8764 (99.9%)                                 ||
## ||             Uniquely mapped : 8756                                         ||
## ||               Multi-mapping : 8                                            ||
## ||                                                                            ||
## ||                    Unmapped : 7                                            ||
## ||                                                                            ||
## ||             Properly paired : 6316                                         ||
## ||         Not properly paired : 2448                                         ||
## ||                   Singleton : 81                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 2                                            ||
## ||  Unexpected fragment length : 2359                                         ||
## ||       Unexpected read order : 6                                            ||
## ||                                                                            ||
## ||                      Indels : 32                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039517_subset_1.fastq                                  ||
## || Input file 2  : SRR1039517_subset_2.fastq                                  ||
## || Output file   : SRR1039517.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:14, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   56% completed, 0.1 mins elapsed, rate=19.4k fragments per second         ||
## ||   Estimated fragment length : 131 bp                                       ||
## ||   80% completed, 0.1 mins elapsed, rate=2.6k fragments per second          ||
## ||   94% completed, 0.1 mins elapsed, rate=2.9k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 11781                                        ||
## ||                      Mapped : 11767 (99.9%)                                ||
## ||             Uniquely mapped : 11760                                        ||
## ||               Multi-mapping : 7                                            ||
## ||                                                                            ||
## ||                    Unmapped : 14                                           ||
## ||                                                                            ||
## ||             Properly paired : 8525                                         ||
## ||         Not properly paired : 3242                                         ||
## ||                   Singleton : 112                                          ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 10                                           ||
## ||  Unexpected fragment length : 3105                                         ||
## ||       Unexpected read order : 15                                           ||
## ||                                                                            ||
## ||                      Indels : 40                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039520_subset_1.fastq                                  ||
## || Input file 2  : SRR1039520_subset_2.fastq                                  ||
## || Output file   : SRR1039520.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:18, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   Estimated fragment length : 143 bp                                       ||
## ||   95% completed, 0.1 mins elapsed, rate=1.6k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 5840                                         ||
## ||                      Mapped : 5829 (99.8%)                                 ||
## ||             Uniquely mapped : 5827                                         ||
## ||               Multi-mapping : 2                                            ||
## ||                                                                            ||
## ||                    Unmapped : 11                                           ||
## ||                                                                            ||
## ||             Properly paired : 4349                                         ||
## ||         Not properly paired : 1480                                         ||
## ||                   Singleton : 50                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 7                                            ||
## ||  Unexpected fragment length : 1418                                         ||
## ||       Unexpected read order : 5                                            ||
## ||                                                                            ||
## ||                      Indels : 36                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //================================= setting ==================================\\
## ||                                                                            ||
## || Function      : Read alignment (RNA-Seq)                                   ||
## || Input file 1  : SRR1039521_subset_1.fastq                                  ||
## || Input file 2  : SRR1039521_subset_2.fastq                                  ||
## || Output file   : SRR1039521.bam (BAM)                                       ||
## || Index name    : homo_sapiens_GRCh38_dna_chromosome_1_rsubread              ||
## ||                                                                            ||
## ||                    ------------------------------------                    ||
## ||                                                                            ||
## ||                               Threads : 1                                  ||
## ||                          Phred offset : 33                                 ||
## ||               # of extracted subreads : 10                                 ||
## ||                        Min read1 vote : 3                                  ||
## ||                        Min read2 vote : 1                                  ||
## ||                     Max fragment size : 600                                ||
## ||                     Min fragment size : 50                                 ||
## ||                        Max mismatches : 3                                  ||
## ||                      Max indel length : 5                                  ||
## ||            Report multi-mapping reads : yes                                ||
## || Max alignments per multi-mapping read : 1                                  ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================ Running (08-Nov-2019 10:16:22, pid=25232) =================\\
## ||                                                                            ||
## || The input file contains base space reads.                                  ||
## || The range of Phred scores observed in the data is [2,41]                   ||
## || Load the 1-th index block...                                               ||
## ||   65% completed, 0.1 mins elapsed, rate=18.5k fragments per second         ||
## ||   Estimated fragment length : 144 bp                                       ||
## ||   83% completed, 0.1 mins elapsed, rate=2.3k fragments per second          ||
## ||   99% completed, 0.1 mins elapsed, rate=2.6k fragments per second          ||
## ||                                                                            ||
## ||                           Completed successfully.                          ||
## ||                                                                            ||
## \\====================================    ====================================//
## 
## //================================   Summary =================================\\
## ||                                                                            ||
## ||             Total fragments : 10163                                        ||
## ||                      Mapped : 10152 (99.9%)                                ||
## ||             Uniquely mapped : 10147                                        ||
## ||               Multi-mapping : 5                                            ||
## ||                                                                            ||
## ||                    Unmapped : 11                                           ||
## ||                                                                            ||
## ||             Properly paired : 7384                                         ||
## ||         Not properly paired : 2768                                         ||
## ||                   Singleton : 87                                           ||
## ||                    Chimeric : 0                                            ||
## ||       Unexpected strandness : 8                                            ||
## ||  Unexpected fragment length : 2661                                         ||
## ||       Unexpected read order : 12                                           ||
## ||                                                                            ||
## ||                      Indels : 41                                           ||
## ||                                                                            ||
## ||                Running time : 0.1 minutes                                  ||
## ||                                                                            ||
## \\============================================================================//
##                                 SRR1039508.bam SRR1039509.bam SRR1039512.bam
## Total_fragments                           7097           7168           8478
## Mapped_fragments                          7090           7157           8461
## Uniquely_mapped_fragments                 7086           7153           8457
## Multi_mapping_fragments                      4              4              4
## Unmapped_fragments                           7             11             17
## Properly_paired_fragments                 5248           5123           6187
## Singleton_fragments                         45             55             58
## More_than_one_chr_fragments                  0              0              0
## Unexpected_strandness_fragments              9             15              8
## Unexpected_template_length                1784           1961           2202
## Inversed_mapping                             4              3              6
## Indels                                      38             28             34
##                                 SRR1039513.bam SRR1039516.bam SRR1039517.bam
## Total_fragments                           7514           8771          11781
## Mapped_fragments                          7507           8764          11767
## Uniquely_mapped_fragments                 7503           8756          11760
## Multi_mapping_fragments                      4              8              7
## Unmapped_fragments                           7              7             14
## Properly_paired_fragments                 5407           6316           8525
## Singleton_fragments                         76             81            112
## More_than_one_chr_fragments                  0              0              0
## Unexpected_strandness_fragments             10              2             10
## Unexpected_template_length                2005           2359           3105
## Inversed_mapping                             9              6             15
## Indels                                      25             32             40
##                                 SRR1039520.bam SRR1039521.bam
## Total_fragments                           5840          10163
## Mapped_fragments                          5829          10152
## Uniquely_mapped_fragments                 5827          10147
## Multi_mapping_fragments                      2              5
## Unmapped_fragments                          11             11
## Properly_paired_fragments                 4349           7384
## Singleton_fragments                         50             87
## More_than_one_chr_fragments                  0              0
## Unexpected_strandness_fragments              7              8
## Unexpected_template_length                1418           2661
## Inversed_mapping                             5             12
## Indels                                      36             41

5 CountTable

fcAirway<-featureCounts(files=bamfls,annot.ext="Homo_sapiens.GRCh38.98.gtf.gz",isGTFAnnotationFile=TRUE, GTF.featureType = "exon", GTF.attrType = "gene_id", useMetaFeatures = TRUE, strandSpecific = 0, 
                    isPairedEnd = TRUE)
## 
##         ==========     _____ _    _ ____  _____  ______          _____  
##         =====         / ____| |  | |  _ \|  __ \|  ____|   /\   |  __ \ 
##           =====      | (___ | |  | | |_) | |__) | |__     /  \  | |  | |
##             ====      \___ \| |  | |  _ <|  _  /|  __|   / /\ \ | |  | |
##               ====    ____) | |__| | |_) | | \ \| |____ / ____ \| |__| |
##         ==========   |_____/ \____/|____/|_|  \_\______/_/    \_\_____/
##        Rsubread 1.34.7
## 
## //========================== featureCounts setting ===========================\\
## ||                                                                            ||
## ||             Input files : 8 BAM files                                      ||
## ||                           o SRR1039508.bam                                 ||
## ||                           o SRR1039509.bam                                 ||
## ||                           o SRR1039512.bam                                 ||
## ||                           o SRR1039513.bam                                 ||
## ||                           o SRR1039516.bam                                 ||
## ||                           o SRR1039517.bam                                 ||
## ||                           o SRR1039520.bam                                 ||
## ||                           o SRR1039521.bam                                 ||
## ||                                                                            ||
## ||              Annotation : Homo_sapiens.GRCh38.98.gtf.gz (GTF)              ||
## ||      Dir for temp files : .                                                ||
## ||                 Threads : 1                                                ||
## ||                   Level : meta-feature level                               ||
## ||              Paired-end : yes                                              ||
## ||      Multimapping reads : counted                                          ||
## || Multi-overlapping reads : not counted                                      ||
## ||   Min overlapping bases : 1                                                ||
## ||                                                                            ||
## ||          Chimeric reads : counted                                          ||
## ||        Both ends mapped : not required                                     ||
## ||                                                                            ||
## \\============================================================================//
## 
## //================================= Running ==================================\\
## ||                                                                            ||
## || Load annotation file Homo_sapiens.GRCh38.98.gtf.gz ...                     ||
## ||    Features : 1371695                                                      ||
## ||    Meta-features : 60623                                                   ||
## ||    Chromosomes/contigs : 47                                                ||
## ||                                                                            ||
## || Process BAM file SRR1039508.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 7097                                                 ||
## ||    Successfully assigned alignments : 5591 (78.8%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039509.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 7168                                                 ||
## ||    Successfully assigned alignments : 6004 (83.8%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039512.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 8478                                                 ||
## ||    Successfully assigned alignments : 6636 (78.3%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039513.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 7514                                                 ||
## ||    Successfully assigned alignments : 6538 (87.0%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039516.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 8771                                                 ||
## ||    Successfully assigned alignments : 6955 (79.3%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039517.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 11781                                                ||
## ||    Successfully assigned alignments : 9804 (83.2%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039520.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 5840                                                 ||
## ||    Successfully assigned alignments : 4581 (78.4%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## || Process BAM file SRR1039521.bam...                                         ||
## ||    Paired-end reads are included.                                          ||
## ||    Total alignments : 10163                                                ||
## ||    Successfully assigned alignments : 8835 (86.9%)                         ||
## ||    Running time : 0.00 minutes                                             ||
## ||                                                                            ||
## ||                                                                            ||
## \\============================================================================//
head(fcAirway$counts)
##                 SRR1039508.bam SRR1039509.bam SRR1039512.bam SRR1039513.bam
## ENSG00000223972              0              0              0              0
## ENSG00000227232              0              0              0              0
## ENSG00000278267              0              0              0              0
## ENSG00000243485              0              0              0              0
## ENSG00000284332              0              0              0              0
## ENSG00000237613              0              0              0              0
##                 SRR1039516.bam SRR1039517.bam SRR1039520.bam SRR1039521.bam
## ENSG00000223972              0              0              0              0
## ENSG00000227232              0              0              0              0
## ENSG00000278267              0              0              0              0
## ENSG00000243485              0              0              0              0
## ENSG00000284332              0              0              0              0
## ENSG00000237613              0              0              0              0
saveRDS(fcAirway,file="fcAirway.rds")
saveRDS(pdata,file="airwayMetaData.rds")