Created r format(Sys.Date(),"%m/%d/%y")
require(immunoSeqR) require(ggplot2) require(gridExtra) require(reshape2) require(knitr) require(rmarkdown) first_run <- TRUE verbose=FALSE r <- function(){render('~/Documents/emj/ImmunoseqResults/immunoSeqR/scripts/rmd/baseline.Rmd', output_file='~/Documents/emj/ImmunoseqResults/new/baseline.pdf', intermediates_dir='~/Documents/emj/ImmunoseqResults/new/')}
ds <- readRDS('~/Documents/emj/ImmunoseqResults/data/baseline/ds_agg.Rds') dict <- readRDS('~/Documents/emj/ImmunoseqResults/data/baseline/dict.Rds') stats <- readRDS('~/Documents/emj/ImmunoseqResults/data/baseline/stats.Rds') load('~/Documents/emj/ImmunoseqResults/data/baseline/olm.Rda') load('~/Documents/emj/ImmunoseqResults/data/baseline/sum.ds.Rda')
\includegraphics[width=\textwidth]{/home/ahopkins/Documents/emj/figures/J0810_timeline.pdf}
The data set contains r sum.ds$nsamp
samples, of r sum.ds$ntypes
types.
There are r format(sum.ds$unique.tcr,big.mark=',')
unique TCRs, collapsed from r
format(sum.ds$total.tcr,big.mark=',')
distinct productive CDR3 sequences. This represents an
average synonymity of r round(sum.ds$avg.syn,3)
. The maximum synonymity (in the parent data
set) was r sum.ds$max.syn
(clone sequence r sum.ds$max.syn.aa
).\par
The available metadata fields in the dictionary are: r paste(names(dict)[-1],collapse=', ')
. The statistics computed are r paste(names(stats)[names(stats)!='fn'],collapse=', ')
.
\clearpage
#put sum of seqs at the end stats <- stats[,c(1,3,2,4)] plot_ds <- merge(dict,stats)
x_val <- 'response' metrics <- names(stats[names(stats)!='fn']) types <- levels(plot_ds$type) for(b in metrics){ l <- vector('list') for(a in types){ l[[a]] <- iseqr_plot_factor(plot_ds,b,x_val,a) } do.call(grid.arrange,c(l,ncol=length(types))) }
\clearpage
metrics <- names(stats[names(stats)!='fn']) l <- vector('list') for(b in metrics){ l[[b]] <- iseqr_plot_factor(plot_ds,b,'type',NA) + theme(axis.text=element_text(size=5)) } do.call(grid.arrange,c(l,ncol=length(types)))
\clearpage
rownames(mm) <- paste0(dict$patient,dict$type) colnames(mm) <- rownames(mm) w <- c(which(dict$response=='NR'),which(dict$response=='R')) tdict <- dict[w,] tmm <- mm[w,w] tmm.m <- melt(tmm) g <- ggplot(tmm.m,aes(x=Var1,y=Var2)) + geom_tile(aes(fill=value)) + coord_fixed() + scale_fill_gradient(low='white',high='black',na.value='white')+ theme(axis.text=element_text(size=5,colour=tdict$response), axis.text.x=element_text(angle=90,hjust=1,vjust=0.5))+ xlab('') + ylab('') + ggtitle('Morisita Overlap') g
\clearpage
#rm(ds) #ds_nt <- readRDS('~/Documents/emj/ImmunoseqResults/mega/fixed/ds_Neoadjuvant_fixed.Rds')
#out <- iseqr_exp_cl(ds_nt,dict) #out <- readRDS('~/Documents/emj/ImmunoseqResults/mega/fixed/exp_cl_out_nadj.Rds')
resp <- sapply(names(out$num_exp),FUN=iseqr_lookup,dict=dict) arm <- sapply(names(out$num_exp),FUN=iseqr_lookup,dict=dict,o_col='arm') df <- data.frame(patient=names(out$num_exp),response=resp,arm=arm,num_exp=out$num_exp) g <- iseqr_plot_factor(df,'num_exp','response',NA) g <- g + ylab('Number of Expanded Clones') f <- iseqr_plot_factor(df,'num_exp','arm',NA) f <- f + ylab('Number of Expanded Clones') grid.arrange(g,f,ncol=2)
\clearpage
mats <- out$total patients <- names(mats) resp <- sapply(patients,FUN=iseqr_lookup,dict=dict) mats <- mats[order(resp)] patients <- patients[order(resp)] resp <- resp[order(resp)] nr <- length(which(resp=='NR')) r <- length(which(resp=='R')) layout <- matrix(nrow=max(r,nr),ncol=2) layout[seq(nr),1] <- seq(nr) layout[,2] <- seq(r)+nr plots <- list() for(a in seq_along(mats)){ tds <- as.data.frame(mats[[a]]) w <- which(dict$patient==patients[a] & dict$type=='PDAC') tds[,4] <- ds_nt[,w] names(tds)[4] <- 'PDAC' #calculate the lfc ttds <- tds ttds[ttds[,1]==0,1] <- 1 ttds[ttds[,2]==0,2] <- 1 tds$lfc <- log(ttds[,2]/ttds[,1],2) rm(ttds) # lowest <- min(tds$p_adj[tds$p_adj>0]) # lowest non-0 p-val tds$p_adj[tds$p_adj==0] <- lowest # clip any p=0 to next lowest tds <- tds[which(tds$PDAC>0),] # restrict to tumor suppressWarnings( l <- max(abs(-log(tds$p_adj,10))) # plot limit ) if(l<1.3){l <- 1.3} # minimum of -log(0.05,10) sh <- shared(tds$PDAC,rowSums(tds[,c(1,2)])) if(nrow(tds)>0){ plots[[a]] <- ggplot(tds,aes(x=-log(p_adj,10),y=PDAC,color=lfc)) + geom_point() + theme_bw()+ ylab('Tumor Counts') + xlab('-log(p)') + ggtitle(paste0(patients[a],'(',resp[a],')',' ',sh,' ','shared')) + scale_colour_gradient2(low='red',mid='black',high='green',midpoint=0)+ scale_x_continuous(limits=c(0,l)) }else{ plots[[a]] <- suppressWarnings(ggplot() + geom_blank()+ ggtitle(paste0(patients[a],'(',resp[a],')',' ',sh,' ','shared')) ) } } g <- grid.arrange(grobs=plots,layout_matrix=layout) #ggsave(g,file='~/Desktop/test_fixed_nadj.pdf',height=2*nrow(layout),width=8)
\clearpage
a <- matrix(c('Arm A:','GVAX','Arm B:','GVAX + IV Cyclophosphamide','Arm C:','GVAX + Oral Cyclophosphamide'),ncol=2,byrow=2) kable(a)
x_val <- 'arm' metrics <- names(stats[names(stats)!='fn']) types <- levels(plot_ds$type) for(b in metrics){ l <- vector('list') for(a in types){ l[[a]] <- iseqr_plot_factor(plot_ds,b,x_val,a) } do.call(grid.arrange,c(l,ncol=length(types))) }
\clearpage
x_val <- 'age' metrics <- names(stats[names(stats)!='fn']) types <- levels(plot_ds$type) for(b in metrics){ l <- vector('list') for(a in types){ l[[a]] <- iseqr_plot_metrics(plot_ds,b,x_val,a,sm=TRUE) } do.call(grid.arrange,c(l,ncol=length(types))) cat(' \n') }
\clearpage
x_val <- 'cells' metrics <- names(stats[names(stats)!='fn']) types <- levels(plot_ds$type) for(b in metrics){ l <- vector('list') for(a in types[types!='PDAC']){ l[[a]] <- iseqr_plot_metrics(plot_ds[!is.na(plot_ds$cells),],b,x_val,a) + theme(axis.text=element_text(size=5)) } do.call(grid.arrange,c(l,ncol=length(types))) cat(' \n') }
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.