In jos4uke/qc4rnaseq: Estimate a posteriori the quality of RNA-seq libraries using count data

\newpage

# treatment of stats file
stats <- paste(opt$stats_dir, basename(opt$stats), sep = .Platform$file.sep)
stats.df <- read.table(stats, header=TRUE, sep="\t")

# treatment of design file
if (!is.null(opt$design)) {
  design <- paste(opt$design_dir, basename(opt$design), sep = .Platform$file.sep)
  design.df <- read.table(design, header=TRUE, sep="\t")
  Factor = design.df[,2]
} 



# dataframe for qplot (library ggplot2)
if (!is.null(opt$design)) {
  datadf <- as.data.frame(cbind(stats.df, Factor))
} else {
  datadf <- as.data.frame(stats.df)
}

Barplot representing the number of raw reads per library

This barplot displays the size of each library (number of raw reads)

Expected results :

• All the libraries should have at least the number of reads warranted by the NGS provider

• The number of reads can be heterogeneous among libraries, due to experimental steps (library preparation, PCR, sequencing ...)

if (!is.null(opt$design))  {
#qplot(data=datadf, x=lib, y=raw_reads.pairs_count, fill= Factor, geom="bar", stat="identity") + theme(axis.text.x = element_text(angle = 90, hjust = 1))
ggplot(data=datadf, aes(x=lib, y=raw_reads.pairs_count)) + geom_bar(stat="identity", aes(fill=Factor)) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
} else {
#qplot(data=datadf, x=lib, y=raw_reads.pairs_count, geom="bar", stat="identity") + theme(axis.text.x = element_text(angle = 90, hjust = 1))
ggplot(data=datadf, aes(x=lib, y=raw_reads.pairs_count)) + geom_bar(stat="identity")  + theme(axis.text.x = element_text(angle = 90, hjust = 1))
}

\newpage

Barplot representing the percentage of mapped_reads / raw_reads (per library)

This barplot displays the % of reads that map on the genome.

Expected results :

• All the libraries should have quite the same % of reads that map on the genome

Application :

• Detect contamination problems (by other organisms)

• Check that the mapping parameters are adequate for all the samples (including the samples that do not come from the reference genome used for the mapping)

if (!is.null(opt$design))  {
#qplot(data=datadf, x=lib, y=100 * mapped_reads.pairs_count/raw_reads.pairs_count, fill= Factor, geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of mapped_reads / raw_reads")
ggplot(data=datadf, aes(x=lib, y=100 * mapped_reads.pairs_count/raw_reads.pairs_count)) + geom_bar(stat="identity", aes(fill=Factor)) + ylim(0,100) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of mapped_reads / raw_reads")
} else {
#qplot(data=datadf, x=lib, y=100 * mapped_reads.pairs_count/raw_reads.pairs_count, geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of mapped_reads / raw_reads")
ggplot(data=datadf, aes(x=lib, y=100 * mapped_reads.pairs_count/raw_reads.pairs_count)) + geom_bar(stat="identity") + ylim(0,100) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of mapped_reads / raw_reads")
}

\newpage

Barplot representing the percentage of specific hits / mapping_hits (per library)

This barplot displays the % of hits that are specific (mapping only once on the genome) among the total number of hits that map on the genome.

Expected results :

• All the libraries should have quite the same % of specific hits

Application :

• Check that the mapping parameters are adequate to ensure specific hits

• Detect contamination problems (for example by ribosomal RNA)

if (!is.null(opt$design))  {
#qplot(data=datadf, x=lib, y=100 * specific_hits/mapping_hits, fill= Factor, geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
ggplot(data=datadf, aes(x=lib, y=100 * specific_hits/mapping_hits)) + geom_bar(stat="identity", aes(fill=Factor)) + ylim(0,100) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
} else {
#qplot(data=datadf, x=lib, y=100 * specific_hits/mapping_hits, geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
ggplot(data=datadf, aes(x=lib, y=100 * specific_hits/mapping_hits)) + geom_bar(stat="identity") + ylim(0,100) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
}

\newpage

Barplot representing the percentage of feature_overlapping_hits / specific_hits (per library)

This barplot displays the % of specific hits that overlap with a feature (gene/RNA).

Expected results :

• All the libraries should have quite the same % of hits that overlap a feature

Application :

• Check the structural annotation

• Check that pair-end libraries are used in the correct orientation

if (!is.null(opt$design))  {
#qplot(data=datadf, x=lib, y=feature_overlapping_hits.specific_hits_percent, fill= Factor, geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of feature_overlapping_hits / specific_hits")
ggplot(data=datadf, aes(x=lib, y=feature_overlapping_hits.specific_hits_percent)) + geom_bar(stat="identity", aes(fill=Factor)) + ylim(0,100) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of feature_overlapping_hits / specific_hits")
} else {
#qplot(data=datadf, x=lib, y=feature_overlapping_hits.specific_hits_percent,  geom="bar", stat="identity", ylim=c(0, 100)) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of feature_overlapping_hits / specific_hits")
ggplot(data=datadf, aes(x=lib, y=feature_overlapping_hits.specific_hits_percent)) + geom_bar(stat="identity") + theme(axis.text.x = element_text(angle = 90, hjust = 1)) + ylab("% of feature_overlapping_hits / specific_hits")
}

jos4uke/qc4rnaseq documentation built on May 19, 2019, 8:45 p.m.