suppressPackageStartupMessages(library(package = "knitr")) suppressPackageStartupMessages(library(package = "Biobase")) suppressPackageStartupMessages(library(package = "tidyverse"))
library (Biobase) library (limma) library (affy) library (stringr) library (dplyr) library (tibble) library (ggplot2) library (pvca) library (ExpressionNormalizationWorkflow) library (sva) library(ggplot2) library(ggthemes) library(reshape2) library(plyr) library (RColorBrewer) library (beeswarm) library (ggbeeswarm) library(gridExtra) library(grid) library(Matching) library(pheatmap) library(janitor) suppressMessages(library(tidyverse))
library (knitr) knitr::opts_chunk$set(echo = TRUE) opts_knit$set (root.dir="~/Dropbox (BCH)/HIPC/IOF/IS2/RDA/CURRENT/CURRENT_2020_05_28")
rm(list=ls()) setwd("~/Dropbox (BCH)/HIPC/IOF/IS2/RDA/CURRENT/CURRENT_2020_05_28") getwd() workDir <- dirname(getwd()) IS2_eset_withResponse_norm_young <- readRDS ("2020_05_28_IS2_eset_withResponse_norm_young.rds") #this one IS2_eset_withResponse_noNorm_young <- readRDS ("2020_05_28_IS2_eset_withResponse_noNorm_young.rds") IS2_eset_noResponse_norm_young <- readRDS ("2020_05_28_IS2_eset_noResponse_norm_young.rds") #this one IS2_eset_noResponse_noNorm_young <- readRDS ("2020_05_28_IS2_eset_noResponse_noNorm_young.rds") IS2_eset_withResponse_norm_older <- readRDS ("2020_05_28_IS2_eset_withResponse_norm_older.rds") #this one IS2_eset_noResponse_norm_older <- readRDS ("2020_05_28_IS2_eset_noResponse_norm_older.rds") #this one IS2_eset_withResponse_noNorm_older <- readRDS ("2020_05_28_IS2_eset_withResponse_noNorm_older.rds") IS2_eset_noResponse_noNorm_older <- readRDS ("2020_05_28_IS2_eset_noResponse_noNorm_older.rds") IS2_immdata <- readRDS ("2020_06_04_IS2_immdata.rds") #this one IS2_immdata_selected <- readRDS ("2020_06_04_IS2_immdata_selected.rds") immdata <- IS2_immdata f<-intersect(featureNames(IS2_eset_noResponse_noNorm_older),featureNames(IS2_eset_noResponse_noNorm_young)) s<-intersect(varLabels(IS2_eset_noResponse_noNorm_older),varLabels(IS2_eset_noResponse_noNorm_young)) IS2_eset<-new("ExpressionSet", exprs=cbind(exprs(IS2_eset_noResponse_noNorm_older[f,]), exprs(IS2_eset_noResponse_noNorm_young[f,])), phenoData=new('AnnotatedDataFrame',rbind(pData(IS2_eset_noResponse_noNorm_older)[,s],pData(IS2_eset_noResponse_noNorm_young)[,s]))) IS2_eset_noResponse_noNorm_young IS2_eset_withResponse_noNorm_young IS2_eset
Figures for Thomas' paper
IS2_all_GE_metaData<- pData (IS2_eset) #Recode IS2_all_GE_metaData$vaccine=recode(IS2_all_GE_metaData$vaccine, "Quadrivalent Meningococcal Vaccine plain-polysaccharide conjugate"="MenACWY-PS") IS2_all_GE_metaData$vaccine=recode(IS2_all_GE_metaData$vaccine, "Quadrivalent Meningococcal Vaccine protein-polysaccharide conjugate"="MenACWY-CRM") IS2_all_GE_metaData$vaccine =recode(IS2_all_GE_metaData$vaccine, "TIV(2011)"="TIV (2011)") IS2_all_GE_metaData$vaccine =recode(IS2_all_GE_metaData$vaccine, "TIV(2013)"="TIV (2013)") IS2_all_GE_metaData$pathogen_vaccinetype=paste0(IS2_all_GE_metaData$pathogen," (",IS2_all_GE_metaData$vaccine_type,")") unique (IS2_all_GE_metaData$pathogen_vaccinetype) #Convert timepoint to numeric IS2_all_GE_metaData$time_post_last_vax=as.numeric(IS2_all_GE_metaData$time_post_last_vax) IS2_all_GE_metaData$time_post_last_vax #Remove saline groups IS2_all_GE_metaData<-IS2_all_GE_metaData%>% filter(vaccine !="Saline") IS2_all_GE_metaData$time_post_last_vax=as.factor(IS2_all_GE_metaData$time_post_last_vax) table (IS2_all_GE_metaData$time_post_last_vax) IS2_all_GE_metaData$time_post_last_vax=revalue(IS2_all_GE_metaData$time_post_last_vax, c("0.166666666666667"="0.167", "21"=">20", "25"=">20", "28"=">20", "37"=">20", "41"=">20", "60"=">20", "70"=">20", "84"=">20")) timepoint <- factor(IS2_all_GE_metaData$time_post_last_vax, levels=sort(as.numeric(unique(IS2_all_GE_metaData$time_post_last_vax)))) #Created variable for the counts to be used for plotting count=rep(1,length(timepoint)) #Combined variables into one dataset IS2_met=as.data.frame(cbind(IS2_all_GE_metaData, timepoint, count))
colourCount= length (unique (IS2_met$pathogen_vaccinetype)) #getPalette = colorRampPalette(brewer.pal(13, "Set3")) getPalette= colorRampPalette (c(brewer.pal(name="Set3", n = 10), brewer.pal(name="Dark2", n = 4))) #ggsave("Fig1BThomasIOF.png", width = 16, height =12, dpi=900) ggplot(IS2_met, aes(x=time_post_last_vax, y=count, fill=pathogen_vaccinetype))+ geom_bar(stat="identity")+ #scale_x_continuous(breaks=seq(0,14))+ #scale_y_continuous(breaks=c(0,250,500,750,1000))+ scale_fill_manual(values=getPalette (colourCount)) + xlab("Day Post-Vaccination")+ #ylab("# of Datasets")+ ylab("Number of Samples")+ labs(fill='Pathogen (Vaccine Type)') + theme_bw() + theme(text = element_text(size=20), axis.text.x = element_text(hjust=1, size=22, angle=45), axis.text.y = element_text(hjust=1, size=22)) + theme (legend.title = element_text(size=14, face="bold"), legend.text = element_text(size = 12)) + theme(legend.key=element_blank(), legend.key.size=unit(10,"point")) + theme (legend.position = "bottom") + # theme(plot.margin=unit(c(1,1,1,1),"cm")) + theme(legend.margin=margin(c(1,2,1,2), unit='cm')) + guides(fill=guide_legend(ncol=2,byrow=TRUE)) + ggtitle ("Summary of Pathogen (Vaccine Type) for Transcriptomics Dataset")+theme(plot.title = element_text(size=16))
Figure 1C
####Focus on Day 0 subjects#### IS2_met$age_reported[is.na(IS2_met$age_reported)] <- 0 IS2_met$age_reported <-as.integer (as.character(IS2_met$age_reported)) table (IS2_met$age_reported) #IS2_all_GE_metaData$age_reported <- cut(IS2_all_GE_metaData$age_reported, breaks=c(-Inf,0, 18, 45, 65, Inf), #labels=c("Not Specified", "Age below 18","Age 18-45","Age 45-65", "Age over 65")) IS2_met$age_imputed[is.na(IS2_met$age_imputed)] <- 0 IS2_met$age_imputed <-as.integer (as.character(IS2_met$age_imputed)) #table (IS2_met$age_imputed) #colnames (IS2_eset) D0<- grepl(pattern = '_0_Days', colnames(IS2_eset)) D0eset <- IS2_eset[, D0] #dim (D0eset) D0pdata <-IS2_met[match(colnames(D0eset), table=IS2_met$uid),] #dim (D0pdata) #[1] 1222 36 #colnames (D0pdata) #table (D0pdata$study_accession) #dim (D0eset) #Features Samples #25860 1222 #options(max.print=1000000) #colnames(D0eset) D0eset <- as.data.frame(D0eset) #Filter those with 0 age for now D0pdata <- D0pdata %>% filter (age_imputed !=0) D0pdata<-D0pdata %>% filter(vaccine != "Saline" ) vaccine <- D0pdata %>% group_by(vaccine) %>% dplyr::summarise(counts = n()) age_reported <-D0pdata$age_reported age_imputed <- D0pdata$age_imputed race <- D0pdata$race gender_imputed <-D0pdata$gender_imputed matrix <- D0pdata$matrix #matrix #head (D0pdata)[1:20] ggplot(D0pdata, aes(x=pathogen_vaccinetype, y=age_imputed)) + geom_boxplot(outlier.shape=NA) + theme_light() + #scale_shape_manual(values=1:nlevels(vaccine)) + theme(axis.text.x = element_text(size=8, angle=60, hjust = 1, face = "bold")) + theme(axis.text.y = element_text(size = 11, hjust = 1, face = "bold")) + scale_y_continuous(breaks = pretty(age_imputed, n = 10)) + ylab("Age") + xlab ("Study")+ #ggtitle("HIPC Signatures 2 Summary: Comparative Meta-analysis Study Metadata") + theme(plot.title = element_text(face="bold", size=11)) + theme(axis.title = element_text(face="bold", size=9)) + scale_fill_manual(values=getPalette (colourCount)) + #scale_color_viridis_d(option="plasma", direction=-1) + geom_jitter(aes(color=vaccine, shape=gender_imputed), position=position_jitter(width=.15, height=.1))+ theme (legend.title = element_text(size=10, face="bold"), legend.text = element_text(size = 10)) + theme(legend.key=element_blank(), legend.key.size=unit(10,"point")) + theme (legend.position = "bottom") #+ # theme(plot.margin=unit(c(1,1,1,1),"cm")) + #theme(legend.margin=margin(c(1,2,1,2), unit='cm'))
Print session info
sessionInfo()
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.