@@ -1,6 +1,6 @@

 Package: limma

-Version: 3.37.10

-Date: 2018-10-15

+Version: 3.37.11

+Date: 2018-10-29

 Title: Linear Models for Microarray Data

 Description: Data analysis, linear models and differential expression for microarray data.

 Author: Gordon Smyth [cre,aut], Yifang Hu [ctb], Matthew Ritchie [ctb], Jeremy Silver [ctb], James Wettenhall [ctb], Davis McCarthy [ctb], Di Wu [ctb], Wei Shi [ctb], Belinda Phipson [ctb], Aaron Lun [ctb], Natalie Thorne [ctb], Alicia Oshlack [ctb], Carolyn de Graaf [ctb], Yunshun Chen [ctb], Mette Langaas [ctb], Egil Ferkingstad [ctb], Marcus Davy [ctb], Francois Pepin [ctb], Dongseok Choi [ctb]

@@ -8,7 +8,7 @@ import(methods)

 importFrom("grDevices", "col2rgb", "dev.cur", "dev.off", "png", "rgb")

 importFrom("graphics", "abline", "axis", "barplot", "coplot", "image",

            "legend", "lines", "matplot", "mtext", "panel.smooth",

-           "par", "plot", "plot.new", "points", "polygon", "rect",

+           "par", "plot", "plot.new", "plot.window", "points", "polygon", "rect",

            "segments", "text", "title")

 importFrom("stats", "approx", "approxfun", "as.dendrogram", "as.dist", "cmdscale",

            "coef", "contr.sum", "contrasts<-", "cov", "cov2cor",

new file mode 100644

@@ -0,0 +1,173 @@

+plotExonJunc <- function(fit, coef=ncol(fit), geneid, genecolname=NULL, FDR=0.05, annotation=NULL)

+#	Assuming the data for diffSplice analysis contains both exons and junctions,

+#	'fit' is a MArrayLM object produced by diffSplice().

+#	To distinguish between exons and junctions, 'fit$genes$Length' are set to 1 for all the junctions.

+#	Since the diffSplice analysis is usually performed after filtering,

+#	the full annotation (e.g. the inbuilt annotation used by featureCounts)

+#	is required for producing the plot.

+#	Yunshun Chen and Gordon Smyth.

+#	Created 9 March 2018. Last modified 29 Oct 2018.

+{

+	if(is.null(genecolname)) 

+		genecolname <- fit$genecolname

+	else

+		genecolname <- as.character(genecolname)

+

+	geneid <- as.character(geneid)

+	i <- fit$genes[, genecolname]==geneid

+	if(!any(i)) stop(paste0(geneid, " not found."))

+

+#	Subsetting

+	fdr <- p.adjust(fit$p.value[,coef], method="BH")

+	fdr <- fdr[i]

+	genes <- fit$genes[i,]

+	strand <- genes$Strand[1]

+	p.value <- fit$p.value[i,coef]

+	coefficients <- fit$coefficients[i,coef]

+

+#	Sorting	exons and junctions by their start positions

+	o <- order(genes$Start, genes$End)

+	genes <- genes[o,]

+	p.value <- p.value[o]

+	coefficients <- coefficients[o]

+	

+#	Which ones are exons? (Junctions are assigned length of 1 prior to the diffSplice analysis)

+	IsExon <- genes$Length > 1L

+	genes.e <- genes[IsExon, ]

+

+#	Check the format of the annotation file.

+	if(!is.null(annotation)){

+		if(is.null(annotation$GeneID)) stop("Annotation file must contain Entrez gene ids.")

+		if(is.null(annotation$Start) | is.null(annotation$End)) stop("Annotation file must contain start-end positions of exons.")

+		if(is.null(annotation$Length)) annotation$Length <- abs(annotation$End - annotation$Start) + 1L

+

+#		Retrieve annotation information for the exons that have been filtered prior to the diffSplice analysis

+		sel <- annotation$GeneID == genes$GeneID[1]

+		genes.e2 <- annotation[sel, ]

+		genes.e2 <- genes.e2[order(genes.e2$Start), ]

+		m <- match(genes.e$Start, genes.e2$Start)

+		genes.e <- genes.e2

+	}

+

+#	Get the start-end positions and the length for all the introns

+	Start.i <- genes.e$End[-nrow(genes.e)] + 1L

+	End.i <- genes.e$Start[-1] - 1L

+	genes.i <- genes.e[-1,]

+	genes.i$Start <- Start.i

+	genes.i$End <- End.i

+	genes.i$Length <- genes.i$End - genes.i$Start + 1L

+

+#	Get the start-end positions for all the junctions

+	if(any(!IsExon)){

+		genes.j <- genes[!IsExon, ]

+		

+#		Extend the plotting range for the gene in case there are junctions outside of the gene range.

+		if(min(genes.j$Start) < min(genes.e$Start)){

+			intron <- genes.i[1,,drop=FALSE]

+			intron$Start <- min(genes.j$Start)

+			intron$End <- min(genes.e$Start) - 1L

+			intron$Length <- intron$End - intron$Start + 1L

+			genes.i <- rbind(intron, genes.i)

+		}

+		if(max(genes.j$End) > max(genes.e$End)){

+			intron <- genes.i[1,,drop=FALSE]

+			intron$Start <- max(genes.e$End) + 1L

+			intron$End <- max(genes.j$End)

+			intron$Length <- intron$End - intron$Start + 1L

+			genes.i <- rbind(genes.i, intron)

+		}

+	}

+

+#	Combine introns and exons for plotting

+	genes.ie <- rbind(cbind(genes.e, Flag="Exon"), cbind(genes.i, Flag="Intron"))

+	genes.ie <- genes.ie[order(genes.ie$Start), ]

+

+#	Scale the length of intron/exon segments for better visualization

+	pseudo.length <- (genes.ie$Length)^.5

+	pseudo.pos <- cumsum((genes.ie$Length)^.5)

+	pseudo.start <- c(0, pseudo.pos[-nrow(genes.ie)])

+	pseudo.end <- pseudo.pos

+	genes.ie <- cbind(genes.ie, pseudo.start=pseudo.start, pseudo.end=pseudo.end, pseudo.length=pseudo.length)

+

+#	Update start-end postions for junctions on the pseudo scale

+	if(any(!IsExon)){

+		genes.j <- cbind(genes.j, pseudo.start=0, pseudo.end=0)

+		for(j in 1:nrow(genes.j)){

+			k <- which(genes.j$Start[j] <= genes.ie$End)[1]

+			genes.j$pseudo.start[j] <- genes.ie$pseudo.end[k] - (genes.ie$End[k] - genes.j$Start[j]) / genes.ie$Length[k] * genes.ie$pseudo.length[k]

+			k <- which(genes.j$End[j] <= genes.ie$End)[1]

+			genes.j$pseudo.end[j] <- genes.ie$pseudo.end[k] - (genes.ie$End[k] - genes.j$End[j]) / genes.ie$Length[k] * genes.ie$pseudo.length[k]

+		}

+	}

+

+#	Setup the plot

+	GeneStart <- min(genes.ie$pseudo.start)

+	GeneEnd <- max(genes.ie$pseudo.end)

+	gene.length <- GeneEnd - GeneStart

+	plot.new()

+	plot.window(xlim=c(GeneStart, GeneEnd), ylim=c(-0.7, 0.7))

+	title(main=paste0(geneid, " (", genes$Strand[1], ")"))

+

+#	Plot gene range

+	rect(xleft=GeneStart, xright=GeneEnd, ybottom=-0.02, ytop=0.02, col="gray", border="gray")

+	if(strand=="+"){

+		tx.left <- "5'"

+		tx.right <- "3'"

+	} else {

+		tx.left <- "3'"

+		tx.right <- "5'"

+	}

+	text(x=-0.02*gene.length, y=0.1, labels=tx.left)

+	text(x=1.02*gene.length, y=0.1, labels=tx.right)

+	

+#	Direction and significance of the diffSplice results

+	up <- coefficients > 0

+	down <- coefficients < 0

+	IsSig <- fdr < FDR

+	#IsSig <- p.adjust(p.value, method="holm") < FDR

+	IsSig.j <- IsSig[!IsExon]

+	down.j <- down[!IsExon]

+	

+#	Colouring

+	col <- rep("black", sum(i))

+	col[up & IsSig] <- "red"

+	col[down & IsSig] <- "dodgerblue"

+	col.e <- col[IsExon]

+	col.j <- col[!IsExon]

+	

+#	Filtered exons are coloured in grey

+	if(!is.null(annotation)){

+		col.e2 <- rep("grey", nrow(genes.e))

+		col.e2[m] <- col.e

+		col.e <- col.e2

+	}

+

+#	Plot exons

+	ex <- genes.ie$Flag=="Exon"

+	rect(xleft=genes.ie$pseudo.start[ex], xright=genes.ie$pseudo.end[ex], ybottom=-0.1,ytop=0.1, col=col.e, border=col.e)

+

+#	Plot junctions

+	if(any(!IsExon)){

+		MidPoint <- (genes.j$pseudo.start + genes.j$pseudo.end)/2

+		y0 <- rep(0.11, sum(!IsExon))

+		y1 <- rep(0.4, sum(!IsExon))

+		y1[IsSig.j] <- 0.6

+		y0[down.j] <- -y0[down.j]

+		y1[down.j] <- -y1[down.j]

+		segments(x0=genes.j$pseudo.start, x1=MidPoint, y0=y0, y1=y1, col=col.j, lwd=2)

+		segments(x0=MidPoint, x1=genes.j$pseudo.end, y0=y1, y1=y0, col=col.j, lwd=2)

+	}

+

+#	Label axis

+	if(genes$Strand[1]=="+")

+		labels <- paste0("Exon.", 1:length(col.e))

+	else

+		labels <- paste0("Exon.", length(col.e):1)

+	axis(side=1, at=(genes.ie$pseudo.start[ex]+genes.ie$pseudo.end[ex])/2, las=2, labels=labels)

+

+	invisible()

+}

+

+

+

+

@@ -3,6 +3,116 @@

 \encoding{UTF-8}

+\section{Version 3.38.0}{\itemize{

+

+\item

+New function plotExonJunc() to plot results from diffSplice().

+

+\item

+New function logsumexp().

+

+\item

+New argument hl.col for volcanoplot(), allowing users to specify

+  the color for the gene names when highlight > 0.

+

+\item

+barcodeplot() no longer assumes that 'statistic'

+  has unique names. Previously it returned an error if

+  names(statistic) contained any duplicated values.

+

+\item

+The colors "blue", "red" and "yellow" used by coolmap() changed to

+  "blue2", "red2" and "yellow2" when used in a color panel with white.

+

+\item

+goana.Rd now explains more explicitly that p-values are unadjusted

+  for multiple testing.

+

+\item

+arrayWeights.Rd now mentions minimum dimensions for expression

+  object.

+

+\item

+More advice on how to choose 'lfc' added to the treat() help page.

+

+\item

+Minor bug fix to the mixed p-value from roast() and mroast() when

+  set.statistic="floormean".

+

+\item

+Bug fix for cumOverlap(), which was under-counting overlaps in some

+  cases.

+

+}}

+

+

+\section{Version 3.36.0}{\itemize{

+

+\item

+New arguments 'quote' and 'row.names' for write.fit(). Write.fit()

+  now outputs row.names by default, which previously were suppressed.

+

+\item

+coolmap() will now accept an arbitrary vector of colors. A new

+  preset "whitered" panel is now also supported.

+

+\item

+In the past, contrasts.fit() always returned NA contrasts for genes

+  with any NA fitted coefficients. contrasts.fit() will now ignore an

+  NA coefficient if the contrast multiplier is zero for that

+  coefficient.

+

+\item

+The treat() default for 'lfc' changed from lfc=0 to lfc=log2(1.2).

+

+\item

+kegga() and goana() now check whether a data.frame has been input

+  by mistake and generates an error. Previously a data.frame value

+  for 'de' was interpreted as a list of gene sets without any error.

+

+\item

+voomWithQualityWeights() now returns 'sample.weights' as a column

+  of the 'targets' data.frame instead of as a separate vector.

+

+\item

+"TestResults" objects now include a 'labels' attribute, defaulting

+  to c("Down","NotSig","Up").

+

+\item

+Functions ebayes() and toptable(), long ago replaced by eBayes() and

+topTable(), are now formally deprecated.

+

+\item

+Update the User's Guide case study that analyses Lrp- E. Coli

+  samples and profiled by Affymetrix arrays.

+

+\item

+Update the Agilent single-channel case study in Section 17.4 of the

+  Users' Guide to use the Agilent gIsWellAboveBG detection column.

+

+\item

+Comments added to voom.Rd and eBayes.Rd to clarify the relationship

+  between limma-trend and voom.

+

+\item

+In the kegga help page, clarify default gene ID system used by

+kegga() when species="dme".

+

+\item

+Update wsva() help page to describe correct number of columns for output.

+

+\item

+bug fix for goana() when 'covariate' is specified and some of the

+  'universe' genes don't have GO annotation.

+

+\item

+Fix to how roast(), mroast() and fry() handle weights. The 'weight' argument

+can now be a matrix, or a vector of length nrow(y), or a vector of length ncol(y).

+This allows the functions to accept precision weights, or gene weights, or sample weights.

+

+}}

+

+

 \section{Version 3.34.0}{\itemize{

 \item 

@@ -62,6 +172,7 @@ plotWithHighlights() checks whether 'status' is a factor and

 Bug fixes for beadCountWeights(). Default is now set correct for

   'design' and the function now works correctly when 'y' is an EList

   object contain bead standard errors but not standard deviations.

+

 }}

@@ -1,3 +1,11 @@

+29 Oct 2018: limma 3.37.11

+

+- New function plotExonJunc() to plot results from diffSplice() when

+  exon-exon junctions as well as exons are included in the count

+  matrix.

+

+- Update NEWS.Rd for Bioconductor 3.8 Release.

+

 14 Oct 2018: limma 3.37.10

 - Edit to barcodeplot() so that it no longer assumes that 'statistic'

new file mode 100644

@@ -0,0 +1,51 @@

+\title{Differential splicing plot with junctions}

+\name{plotExonJunc}

+\alias{plotExonJunc}

+\description{

+Plot differential usage results by exons and junctions for the specified gene and highlight the significantly spliced exons and junctions.

+}

+\usage{

+plotExonJunc(fit, coef=ncol(fit), geneid, genecolname=NULL, FDR=0.05, annotation=NULL)

+}

+\arguments{

+  \item{fit}{\code{MArrayLM} fit object produced by \code{diffSplice}. Must have the Entrez gene ids for all the exons and junctions stored in \code{fit$genes$GeneID}, length information for all the exons and junctions stored in \code{fit$genes$Length} and the strand information stored in \code{fit$genes$Strand}. To distinguish between exons and junctions, \code{fit$genes$Length} has to be set to 1 for all the junctions.}

+  \item{coef}{the coefficient (column) of fit for which differentially splicing is assessed.}

+  \item{geneid}{character string, ID of the gene to plot.}

+  \item{genecolname}{column name of \code{fit$genes} containing \code{geneid}.}

+  \item{FDR}{numeric, highlight exons and junctions with false discovery rate less than this cutoff. Red indicates up-regulation whereas blue indicates down-regulation. The FDR of the individual exon/junction is calculated based on the exon-level t-statistics test for differences between each exon/junction and all other exons/junctions for the same gene.}

+  \item{annotation}{data frame containing the full exon annotation of the corresponding species. Must have the Entrez gene ids for all the exons stored in the \code{GeneID} column, start and end positions for all the exons stored in the \code{Start} and \code{End} columns, respectively.}

+}

+

+\details{

+Plot differential usage results by exons and junctions for the specified gene.

+The significantly spliced individual exons are highlighted as red blocks if up-regulated and blue blocks if down-regulated. 

+All other exons are displayed as black blocks.

+The significantly spliced individual junctions are highlighted as red lines if up-regulated and blue lines if down-regulated.

+All other junctions are displayed as black lines.

+

+Since the \code{diffSplice} analysis is usually performed after filtering, the full annotation (e.g. the inbuilt annotation in \code{featureCounts}) is highly recommended for producing the plot. When \code{annotation} is provided, the filtered exons are displayed as grey blocks.

+}

+

+\value{A plot is created on the current graphics device.}

+\author{Yunshun Chen and Gordon Smyth}

+\seealso{

+\code{\link{diffSplice}}, \code{\link{topSplice}}

+}

+\examples{

+\dontrun{

+# diffSplice analysis

+v <- voom(dge, design)

+fit <- lmFit(v, design)

+ex <- diffSplice(fit, geneid="GeneID")

+

+# Get full annotation from Rsubread

+library(Rsubread)

+annotation.full <- getInBuiltAnnotation("mm10")

+

+# Make a plot

+plotExonJunc(ex, geneid="Foxp1", genecolname="Symbol", annotation=annotation.full)

+}

+}

+

+\keyword{hplot}

+\keyword{rna-seq}