Browse code

Commit made by the Bioconductor Git-SVN bridge.

Commit id: 48b838ee337e94cff1f69151d75091b17bd7257f

msaplot


Commit id: 6083e201a71a1ff5ae2d11e635fb91e9e9c0b3fd

update appveyor.yml



git-svn-id: https://hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/ggtree@104099 bc3139a8-67e5-0310-9ffc-ced21a209358

g.yu authored on 22/05/2015 08:47:09
Showing 7 changed files

... ...
@@ -1,7 +1,7 @@
1 1
 Package: ggtree
2 2
 Type: Package
3 3
 Title: a phylogenetic tree viewer for different types of tree annotations
4
-Version: 1.1.6
4
+Version: 1.1.7
5 5
 Author: Guangchuang Yu and Tommy Tsan-Yuk Lam
6 6
 Maintainer: Guangchuang Yu <guangchuangyu@gmail.com>
7 7
 Description: ggtree extends the ggplot2 plotting system which implemented the
... ...
@@ -48,6 +48,7 @@ export(groupOTU)
48 48
 export(gzoom)
49 49
 export(hilight)
50 50
 export(merge_tree)
51
+export(msaplot)
51 52
 export(phylopic)
52 53
 export(plot)
53 54
 export(read.baseml)
... ...
@@ -114,6 +115,7 @@ importFrom(ggplot2,element_rect)
114 115
 importFrom(ggplot2,element_text)
115 116
 importFrom(ggplot2,fortify)
116 117
 importFrom(ggplot2,geom_point)
118
+importFrom(ggplot2,geom_rect)
117 119
 importFrom(ggplot2,geom_segment)
118 120
 importFrom(ggplot2,geom_text)
119 121
 importFrom(ggplot2,geom_tile)
... ...
@@ -124,6 +126,7 @@ importFrom(ggplot2,labs)
124 126
 importFrom(ggplot2,scale_color_manual)
125 127
 importFrom(ggplot2,scale_fill_discrete)
126 128
 importFrom(ggplot2,scale_fill_gradient)
129
+importFrom(ggplot2,scale_fill_manual)
127 130
 importFrom(ggplot2,scale_x_continuous)
128 131
 importFrom(ggplot2,scale_x_reverse)
129 132
 importFrom(ggplot2,scale_y_continuous)
... ...
@@ -1,3 +1,7 @@
1
+CHANGES IN VERSION 1.1.7
2
+------------------------
3
+ o msaplot for adding multiple sequence alignment <2015-05-22, Fri>
4
+ 
1 5
 CHANGES IN VERSION 1.1.6
2 6
 ------------------------
3 7
  o add vertical_only parameter to scaleClade and set to TRUE by default.
... ...
@@ -61,6 +61,92 @@ gheatmap <- function(p, data, offset=0, width=NULL, low="green", high="red",
61 61
     return(p2)
62 62
 }
63 63
 
64
+##' multiple sequence alignment with phylogenetic tree
65
+##'
66
+##' 
67
+##' @title msaplot
68
+##' @param p tree view
69
+##' @param fasta fasta file, multiple sequence alignment
70
+##' @param offset offset of MSA to tree
71
+##' @param width width of each character
72
+##' @param color color 
73
+##' @param window specific a slice to display
74
+##' @return tree view
75
+##' @export
76
+##' @importFrom Biostrings readBStringSet
77
+##' @importFrom colorspace rainbow_hcl
78
+##' @importFrom ggplot2 geom_segment
79
+##' @importFrom ggplot2 geom_rect
80
+##' @importFrom ggplot2 scale_fill_manual
81
+##' @author Guangchuang Yu
82
+msaplot <- function(p, fasta, offset=0, width=NULL, color=NULL, window=NULL){
83
+    aln <- readBStringSet(fasta)
84
+    if (is.null(window)) {
85
+        window <- c(1, width(aln)[1])
86
+    }
87
+    slice <- seq(window[1], window[2], by=1)
88
+    
89
+    seqs <- lapply(1:length(aln), function(i) {
90
+        x <- toString(aln[i])
91
+        seq <- substring(x, slice, slice)
92
+
93
+        seq[seq == '?'] <- '-'
94
+        seq[seq == '*'] <- '-'
95
+        seq[seq == ' '] <- '-'
96
+        return(seq)
97
+    })
98
+    names(seqs) <- names(aln)
99
+    
100
+    if(is.null(color)) {
101
+        alphabet <- unlist(seqs) %>% unique
102
+        alphabet <- alphabet[alphabet != '-']
103
+        color <- rainbow_hcl(length(alphabet))
104
+        names(color) <- alphabet
105
+        color <- c(color, '-'=NA)
106
+    }
107
+
108
+    df <- p$data
109
+    if (is.null(width)) {
110
+        width <- (df$x %>% range %>% diff)/500
111
+    }
112
+
113
+    df=df[df$isTip,]
114
+    start <- max(df$x) * 1.02 + offset
115
+
116
+    seqs <- seqs[df$label[order(df$y)]]
117
+    ## seqs.df <- do.call("rbind", seqs)
118
+
119
+    h <- ceiling(diff(range(df$y))/length(df$y))
120
+    xmax <- start + seq_along(slice) * width
121
+    xmin <- xmax -width
122
+    y <- sort(df$y)
123
+    ymin <- y - 0.4 *h
124
+    ymax <- y + 0.4 *h
125
+
126
+    from <- to <- NULL
127
+    
128
+    lines.df <- data.frame(from=min(xmin), to=max(xmax), y = y)
129
+
130
+    p <- p + geom_segment(data=lines.df, aes(x=from, xend=to, y=y, yend=y))
131
+    msa <- lapply(1:length(y), function(i) {
132
+        data.frame(name=names(seqs)[i],
133
+                   xmin=xmin,
134
+                   xmax=xmax,
135
+                   ymin=ymin[i],
136
+                   ymax=ymax[i],
137
+                   seq=seqs[[i]])
138
+    })
139
+
140
+    msa.df <- do.call("rbind", msa)
141
+
142
+    p <- p + geom_rect(data=msa.df, aes(x=xmin, y=ymin, 
143
+                           xmin=xmin, xmax=xmax,
144
+                           ymin=ymin, ymax=ymax, fill=seq)) +
145
+                               scale_fill_manual(values=color)
146
+
147
+    return(p)
148
+}
149
+
64 150
 ##' scale x for tree with heatmap
65 151
 ##'
66 152
 ##' 
... ...
@@ -6,7 +6,7 @@
6 6
   and compare the feature from parent to children and add store the info in tree object.
7 7
   So with ancestral sequences inferred by HYPHY or PAML, we can support any type of sequence feature comparison,
8 8
   not only substitution supported internally.
9
-+ support more features that can be plotted at the right hand side of the tree,
10
-  not only matrix but also for example sequence alignment or specific window of alignment.
9
++ support more features that can be plotted at the right hand side of the tree.
10
+  - ~~multiple sequence alignment~~ _Now implemented in 1.1.7_
11 11
   
12 12
 
... ...
@@ -13,9 +13,7 @@ install:
13 13
 
14 14
 build_script:
15 15
   - travis-tool.sh install_deps
16
-  - travis-tool.sh install_bioc EBImage
17
-  - travis-tool.sh install_bioc BiocStyle
18
-  - travis-tool.sh install_bioc Biostrings
16
+  - travis-tool.sh install_bioc BiocStyle Biostrings EBImage
19 17
 ##  - travis-tool.sh install_bioc_deps
20 18
   
21 19
 test_script:
22 20
new file mode 100644
... ...
@@ -0,0 +1,31 @@
1
+% Generated by roxygen2 (4.1.1): do not edit by hand
2
+% Please edit documentation in R/gplot.R
3
+\name{msaplot}
4
+\alias{msaplot}
5
+\title{msaplot}
6
+\usage{
7
+msaplot(p, fasta, offset = 0, width = NULL, color = NULL, window = NULL)
8
+}
9
+\arguments{
10
+\item{p}{tree view}
11
+
12
+\item{fasta}{fasta file, multiple sequence alignment}
13
+
14
+\item{offset}{offset of MSA to tree}
15
+
16
+\item{width}{width of each character}
17
+
18
+\item{color}{color}
19
+
20
+\item{window}{specific a slice to display}
21
+}
22
+\value{
23
+tree view
24
+}
25
+\description{
26
+multiple sequence alignment with phylogenetic tree
27
+}
28
+\author{
29
+Guangchuang Yu
30
+}
31
+