In ahopki14/immunoSeqR: Tools for analyzing data from Adaptive's ImmunoSeq platform

Neoadjuvant Study

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/neoadjuvant.Rmd',
output_file='~/Documents/emj/ImmunoseqResults/new/neoadjuvant.pdf',
intermediates_dir='~/Documents/emj/ImmunoseqResults/new/')}

dict <- readRDS('~/Documents/emj/ImmunoseqResults/data/neoadjuvant/dict.Rds')
stats <- readRDS('~/Documents/emj/ImmunoseqResults/data/neoadjuvant/stats.Rds')
plot_ds <- readRDS('~/Documents/emj/ImmunoseqResults/data/neoadjuvant/plot_ds.Rds')
load('~/Documents/emj/ImmunoseqResults/data/neoadjuvant/olm.Rda')
load('~/Documents/emj/ImmunoseqResults/data/neoadjuvant/sum.ds.Rda')

Background

\includegraphics[width=\textwidth]{/home/ahopkins/Documents/emj/figures/J0810_timeline.pdf}

Dataset

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

Metadata

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

By Response

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

Change in Metrics

w <- grep('Log',names(plot_ds))
l <- vector('list')
for(a in seq(length(w))){
  l[[a]] <- iseqr_plot_factor(plot_ds,names(plot_ds)[w[a]],'response','POST') +
                              geom_hline(yintercept=0,alpha=0.1) + ggtitle('PRE to POST')
}
do.call(grid.arrange,c(l,ncol=length(w)))
````
\clearpage 

### By Type
```r
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

Overlaps

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

Expanded Clones

w <- grep('Number of Expanded',names(plot_ds))
g <- iseqr_plot_factor(plot_ds,names(plot_ds)[w],'response',type='POST') 
f <- iseqr_plot_factor(plot_ds,names(plot_ds)[w],'arm',type='POST')
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

Confounders

Treatment Arm

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

Age

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

Cell Number

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')
}

ahopki14/immunoSeqR documentation built on May 7, 2019, 2:54 a.m.