Я пытаюсь настроить R shiny app, который позволит просматривать три типа графиков, относящихся к данным экспрессии генов.
Данные состоят из:
A data.frame, который имеет выходные данные анализа дифференциальной экспрессии (каждая строка представляет собой ген, а столбцы представляют размеры эффекта и их p-значения):
set.seed(1)
model.df <- data.frame(id = paste0("g",1:30),symbol = sample(LETTERS[1:5],30,replace=T),
group.effect.size = rnorm(30), group.p.value = runif(30,0,1),
sex.effect.size = rnorm(30), sex.p.value = runif(30,0,1),
stringsAsFactors = F)
A data.frame, который имеетплан исследования (каждая строка представляет собой образец, а столбцы - это факторы, с которыми связана выборка):
set.seed(1)
design.df <- data.frame(group = c(rep("A",6),rep("B",6)), sex = rep(c(rep("F",3),rep("M",3)),2), replicate = rep(1:6,2)) %>%
dplyr::mutate(sample = paste0(group,".",sex,"_",replicate))
design.df$group <- factor(design.df$group, levels = c("A","B"))
design.df$sex <- factor(design.df$sex, levels = c("F","M"))
A matrix, который имеет численность (каждая строка представляет собой ген, а каждая -столбец является образцом):
set.seed(1)
abundance.mat <- matrix(rnorm(30*12), nrow=30, ncol=12, dimnames=list(model.df$id,design.df$sample))
A data.frame, который имеет результаты анализа обогащения набора генов (каждая строка является именем набора, а столбцы являются тестом на обогащение p-значения для каждого фактора в design.df):
set.seed(1)
gsea.df <- data.frame(set.name = paste0("S",1:4), group.p.value = format(round(runif(4,0,1),2),scientific = T), sex.p.value = format(round(runif(4,0,1),2),scientific = T), stringsAsFactors = F)
И, наконец, data.frame, который связывает гены с каждым set.name в gsea.df:
set.seed(1)
gene.sets.df <- do.call(rbind,lapply(1:4,function(s) data.frame(set.name = paste0("S",s), id = sample(model.df$id,10,replace = F),stringsAsFactors = F)))
Я бы хотел, чтобы shiny app включил просмотр следующих типов графиков:
Feature Plot - сюжетУровень экспрессии одного выбранного пользователем гена на y-axis и образца на x-axis, который будет объединен со вставкой гусеничного графика, показывающего предполагаемые эффекты: 
Feature User-Defined Sets Plot - то же самое, что и Feature Plot, но вместо того, чтобы показывать один-единственный выбранный ген, он покажет набор выбранных пользователем генов и, следовательно, вместо точек он покажет скрипки распределений: 
Feature Sets GSEA Plot - объединенный список графиков вулканов, где на каждом из них x-axis - это величина эффекта фактора, y-axis - -log10 (p-значение) эффекта, и гены окрашены в красный цвет, если они принадлежат выбранному набору генов: 
Вот три функции для генерации этих фигур:выбор пользователя:
featurePlot <- function(selected.id)
{
replicate.df <- reshape2::melt(abundance.mat[which(rownames(abundance.mat) == selected.id),,drop=F], varnames=c("id","sample")) %>%
dplyr::left_join(design.df)
effects.df <- data.frame(factor.name = c("group","sex"),
effect.size = c(dplyr::filter(model.df,id == selected.id)$group.effect.size,dplyr::filter(model.df,id == selected.id)$sex.effect.size),
p.value = c(dplyr::filter(model.df,id == selected.id)$group.p.value,dplyr::filter(model.df,id == selected.id)$sex.p.value),
stringsAsFactors = F)
effects.df$factor.name <- factor(effects.df$factor.name, levels = c("group","sex"))
main.plot <- ggplot(replicate.df,aes(x=replicate,y=value,color=group,shape=sex))+
geom_point(size=3)+facet_grid(~group,scales="free_x")+
labs(x="Replicate",y="TPM")+theme_minimal()
xlims <- c(-1*max(abs(effects.df$effect.size))-0.1*max(abs(effects.df$effect.size)),max(abs(effects.df$effect.size))+0.1*max(abs(effects.df$effect.size)))
effects.plot <- ggplot(effects.df,aes(x=effect.size,y=factor.name,color=factor.name))+
geom_point()+
geom_vline(xintercept=0,linetype="longdash",colour="black",size=0.25)+theme_minimal()+xlim(xlims)+
theme(legend.position="none")+ylab("")+xlab("Effect Size")
null.plot <- ggplot(data.frame())+geom_point()+geom_blank()+theme_minimal()
combined.plot <- gridExtra::arrangeGrob(main.plot,gridExtra::arrangeGrob(null.plot,effects.plot,ncol=1),nrow=1,ncol=2,widths=c(5,2.5))
return(combined.plot)
}
featureSetPlot <- function(selected.ids)
{
replicate.df <- reshape2::melt(abundance.mat[which(rownames(abundance.mat) %in% selected.ids),,drop=F], varnames=c("id","sample")) %>%
dplyr::left_join(design.df)
replicate.df$replicate <- as.factor(replicate.df$replicate)
effects.df <- data.frame(factor.name = c("group","sex"),
effect.size = c(dplyr::filter(model.df,id %in% selected.ids)$group.effect.size,dplyr::filter(model.df,id %in% selected.ids)$sex.effect.size),
p.value = c(dplyr::filter(model.df,id %in% selected.ids)$group.p.value,dplyr::filter(model.df,id %in% selected.ids)$sex.p.value),
stringsAsFactors = F)
effects.df$factor.name <- factor(effects.df$factor.name, levels = c("group","sex"))
main.plot <- ggplot(replicate.df,aes(x=replicate,y=value,color=group,fill=sex))+
geom_violin(trim=F,draw_quantiles=c(0.25,0.5,0.75),alpha=0.25)+facet_grid(~group,scales="free_x")+
labs(x="Replicate",y="TPM")+theme_minimal()
effects.plot <- ggplot(effects.df,aes(y=effect.size,x=factor.name,color=factor.name,fill=factor.name))+
geom_violin(trim=F,draw_quantiles=c(0.25,0.5,0.75),alpha=0.25)+coord_flip()+
geom_hline(yintercept=0,linetype="longdash",colour="black",size=0.25)+theme_minimal()+
theme(legend.position="none")+xlab("")+ylab("Effect Size Distribution")
null.plot <- ggplot(data.frame())+geom_point()+geom_blank()+theme_minimal()
combined.plot <- gridExtra::arrangeGrob(main.plot,gridExtra::arrangeGrob(null.plot,effects.plot,ncol=1),nrow=1,ncol=2,widths=c(5,2.5))
return(combined.plot)
}
gseaPlot <- function(selected.set)
{
plot.df <- model.df %>%
dplyr::left_join(gene.sets.df %>% dplyr::filter(set.name == selected.set))
plot.df$set.name[which(is.na(plot.df$set.name))] <- "non.selected"
plot.df$set.name <- factor(plot.df$set.name, levels = c("non.selected",selected.set))
factor.names <- c("group","sex")
gsea.volcano.plot <- lapply(factor.names,function(f)
plotly::plot_ly(type='scatter',mode="markers",marker=list(size=5),color=plot.df$set.name,colors=c("lightgray","darkred"),x=plot.df[,paste0(f,".effect.size")],y=-log10(plot.df[,paste0(f,".p.value")]),showlegend=F) %>%
plotly::layout(annotations=list(showarrow=F,x=0.5,y=0.95,align="center",xref="paper",xanchor="center",yref="paper",yanchor="bottom",font=list(size=12,color="darkred"),text=paste0(f," (",dplyr::filter(gsea.df,set.name == selected.set)[,paste0(f,".p.value")],")")),
xaxis=list(title=paste0(f," Effect"),zeroline=F),yaxis=list(title="-log10(p-value)",zeroline=F))
) %>% plotly::subplot(nrows=1,shareX=F,shareY=T,titleX=T,titleY=T) %>%
plotly::layout(title=selected.set)
return(gsea.volcano.plot)
}
Таким образом:
plot.type.choices <- c('Feature User-Defined Set Plot','Feature Sets GSEA Plot','Feature Plot')
Итак, первые два используют ggplot2 для генерации каждой из двух фигур, которые они объединяют, что затем достигается с помощью gridExtra::arrangeGrob,Последний использует plotly.
Вот часть кода shiny, которую я пробовал, но без удачи:
server <- function(input, output)
{
out.plot <- reactive({
if(input$plotType == "Feature Plot"){
out.plot <- featurePlot(selected.id=dplyr::filter(model.df,symbol == input$symbol)$id[1])
} else if(input$plotType == "Feature User-Defined Set Plot"){
out.plot <- featureSetPlot(selected.ids=unique(dplyr::filter(model.df,symbol == input$set.symbols)$id))
} else if(input$plotType == "Feature Sets GSEA Plot"){
out.plot <- gseaVolcanoPlot(selected.set=input$set.name)
}
})
output$out.plot <- renderPlot({
if(input$plotType != "Feature Sets GSEA Plot"){
grid::grid.draw(out.plot())
} else{
out.plot()
}
})
output$save <- downloadHandler(
filename = function() {
paste0("./plot.pdf")
},
content = function(file) {
ggsave(out.plot(),filename=file,width=10,height=5)
}
)
}
ui <- fluidPage(
tags$style(type="text/css",".shiny-output-error { visibility: hidden; }",".shiny-output-error:before { visibility: hidden; }"),
titlePanel("Results Explorer"),
sidebarLayout(
sidebarPanel(
# select plot type
selectInput("plotType","Plot Type",choices=plot.type.choices),
#in case Feature User-Defined Set Plot was chosen select the genes
conditionalPanel(condition="input.plotType=='Feature User-Defined Set Plot'",
selectizeInput(inputId="set.symbols",label="Features Set Symbols",choices=unique(model.df$symbol),selected=model.df$symbol[1],multiple=T)),
#in case Feature Sets GSEA Plot was chosen select the databses
conditionalPanel(condition="input.plotType=='Feature Sets GSEA Plot'",
selectizeInput(inputId="set.name",label="Set Name",choices=unique(gene.sets.df$set.name),selected=gene.sets.df$set.name[1],multiple=F)),
#in case Feature Plot was chosen select the gene
conditionalPanel(condition="input.plotType=='Feature Plot'",
selectizeInput(inputId="symbol",label="Feature Symbol",choices=unique(model.df$symbol),selected=unique(model.df$symbol)[1],multiple=F)),
downloadButton('save', 'Save to File')
),
mainPanel(
plotOutput("output.plot")
)
)
)
shinyApp(ui = ui, server = server)
Я подозреваю, что renderPlot здесь может быть проблема, поскольку мне, вероятно, придется использовать plotly::renderPlotly для опции Feature Sets GSEA Plot, но я не совсем уверен, как связать все это в части shiny server.
Другое существующее осложнение, и было бы неплохо найти решение, заключается в том, что генные символы не являются уникальными идентификаторами гена WRT (как показано в model.df).Поэтому было бы неплохо добавить список, который добавится, если пользователь выберет опцию Feature Plot, и в этом списке будет показано подмножество идентификаторов генов, которым соответствует выбранный символ (dplyr::filter(model.df == input$symbol)$id)
Спасибо!