Commit 00a02591 authored by Hanna  JULIENNE's avatar Hanna JULIENNE
Browse files

add plot for all simulation tested

parent f083000c
......@@ -149,4 +149,4 @@ if __name__ == '__main__':
Ssize.N_effective = Ssize.N_effective.astype(int, copy=False)
Ssize.to_csv("../../data/external/meta_data/N_effective.csv")
print(Ssize.head())
simulate_Z_scores("../../data/processed/Simulated/Genotypes/genotype2.csv", "../../data/processed/Simulated/Genotypes/LD_matrix2.csv", Ssize,"../../data/processed/Simulated/Zscores/Zscore", amp=0.05)
simulate_Z_scores("../../data/processed/Simulated/Genotypes/genotype2.csv", "../../data/processed/Simulated/Genotypes/LD_matrix2.csv", Ssize,"../../data/processed/Simulated/Zscores/Zscore", amp=0.25)
......@@ -5,38 +5,39 @@ library(cowplot)
setwd("/mnt/zeus/GGS/PROJECT_imputation_covidhg/hgcovid_imputation/src/visualization")
N_eff = fread("../../data/external/meta_data/N_effective.csv")
for(cohort in c("", "_small_cohort")){
N_eff = fread(paste0("../../data/external/meta_data/N_effective",cohort,".csv"))
N_eff = as.data.frame(N_eff)
row.names(N_eff) = paste0("Z_",N_eff$study)
print(head(N_eff))
Imputation = fread("../../data/processed/Simulated/Imputed/Imputed_one_causal.csv")
for( tag in c("null","one_causal", 'two_opposite', 'two_causal')){
imp_file = paste0("../../data/processed/Simulated/Imputed/Imputed",cohort,"_",tag,".csv")
Imputation = fread(imp_file)
Imputation
cor_ref = data.frame(correlation = cor(Imputation[,which(grepl("Z_",names(Imputation))), with=FALSE])[,1])
cor_ref["N"] = N_eff[row.names(cor_ref), "N_effective"]
head(cor_ref)
p_cor = ggplot(cor_ref, aes(x=N, y=correlation)) + geom_line() + geom_point() + theme_minimal()
N_eff = as.data.frame(N_eff)
row.names(N_eff) = paste0("Z_",N_eff$Name)
N_eff[order(N_eff$N_effective),]
cor_ref = data.frame(correlation = cor(Imputation[,which(grepl("Z_",names(Imputation))), with=FALSE])[,1])
cor_ref
cor_ref["N"] = N_eff[row.names(cor_ref), "N_effective"]
beta_hat = melt(Imputation[,c(1,which(grepl("Beta_",names(Imputation)))), with=FALSE], id.vars=1)
p_cor = ggplot(cor_ref, aes(x=N, y=correlation)) + geom_line() + geom_point() + theme_minimal()
p_cor
beta_hat = melt(Imputation[,c(1,which(grepl("Beta_",names(Imputation)))), with=FALSE], id.vars=1)
beta_hat
psig = ggplot(beta_hat, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat[variable=="Beta_ref",], color="midnightblue",lwd=1.1)
psig =psig + theme_minimal() + xlab("snp") + ylab("beta")
psig
beta_hat_prec = melt(Imputation[,c("V1","Beta_ref",paste0("Beta_",N_eff[N_eff$N_effective > 100,"Name"])), with=FALSE], id.vars=1)
psig_prec = ggplot(beta_hat_prec, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat_prec[variable=="Beta_ref",], color="midnightblue",lwd=1.1)
psig_prec = psig_prec + theme_minimal() + xlab("snp") + ylab("beta")
psig = ggplot(beta_hat, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat[variable=="Beta_ref",], color="midnightblue",lwd=1.1)
psig =psig + theme_minimal() + xlab("snp") + ylab("beta")
psig_prec
Beta_scatter = melt(Imputation[,c(1,which(grepl("Beta_",names(Imputation)))), with=FALSE], id.vars=c(1,2))
ggplot(Beta_scatter, aes(x=Beta_ref, y=value))+geom_point()
beta_hat_prec = melt(Imputation[,c("V1","Beta_ref",paste0("Beta_",N_eff[N_eff$N_effective > 50,"study"])), with=FALSE], id.vars=1)
psig_prec = ggplot(beta_hat_prec, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat_prec[variable=="Beta_ref",], color="midnightblue",lwd=1.1)
psig_prec = psig_prec + theme_minimal() + xlab("snp") + ylab("beta")
Beta_scatter = melt(Imputation[,c(1,which(grepl("Beta_",names(Imputation)))), with=FALSE], id.vars=c(1,2))
ggplot(Beta_scatter, aes(x=Beta_ref, y=value))+geom_point()
Beta_scatter_prec = melt(Imputation[,c("V1","Beta_ref",paste0("Beta_",N_eff[N_eff$N_effective > 100,"Name"])), with=FALSE], id.vars=c(1,2))
Beta_scatter_prec = ggplot(Beta_scatter_prec, aes(x=Beta_ref, y=value, color=variable))+geom_point() + scale_colour_hue() + theme_minimal()+theme(legend.pos="none") + xlab('Beta') + ylab("imputed Beta")
png("../../reports/figures/Signal_variability_imputed.png", width=900, height=325)
plot_grid(p_cor, psig_prec, Beta_scatter_prec, labels=c("A", "B", "C"), nrow=1)
dev.off()
Beta_scatter_prec = melt(Imputation[,c("V1","Beta_ref",paste0("Beta_",N_eff[N_eff$N_effective > 50,"study"])), with=FALSE], id.vars=c(1,2))
Beta_scatter_prec = ggplot(Beta_scatter_prec, aes(x=Beta_ref, y=value, color=variable))+geom_point() + scale_colour_hue() + theme_minimal()+theme(legend.pos="none") + xlab('Beta') + ylab("imputed Beta")
panel = plot_grid(p_cor, psig_prec, Beta_scatter_prec, labels=c("A", "B", "C"), nrow=1)
save_plot(paste0("../../reports/figures/Signal_variability_imputed",cohort,"_",tag,".png"), panel, base_asp=1.8, bg='white')
}
}
......@@ -3,27 +3,47 @@ library(ggplot2)
library(cowplot)
setwd("/mnt/zeus/GGS/PROJECT_imputation_covidhg/hgcovid_imputation/src/visualization")
setwd("/pasteur/zeus/projets/p02/GGS_WKD/PROJECT_imputation_covidhg/hgcovid_imputation/src/visualization")
N_eff = fread("../../data/external/meta_data/N_effective.csv")
Zscores = fread("../../data/processed/Simulated/Zscores/Zscore_one_causal.csv")
for(cohort in c("", "_small_cohort")){
N_eff = fread(paste0("../../data/external/meta_data/N_effective",cohort,".csv"))
N_eff = as.data.frame(N_eff)
N_eff = as.data.frame(N_eff)
row.names(N_eff) = paste0("Z_",N_eff$Name)
row.names(N_eff) = paste0("Z_",N_eff$study)
print(head(N_eff))
cor_ref = data.frame(correlation = cor(Zscores[,which(grepl("Z_",names(Zscores))), with=FALSE])[,1])
cor_ref["N"] = N_eff[row.names(cor_ref), "N_effective"]
for( tag in c("null","one_causal", 'two_opposite', 'two_causal')){
print(cohort)
print(tag)
Zscores = fread(paste0("../../data/processed/Simulated/Zscores/Zscore",cohort,"_",tag,".csv"))
p_cor = ggplot(cor_ref, aes(x=N, y=correlation)) + geom_line() + geom_point() + theme_minimal()
beta_hat = melt(Zscores[,c(1,which(grepl("Beta_",names(Zscores)))), with=FALSE], id.vars=1)
psig = ggplot(beta_hat, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat[variable=="Beta_meta-analysis",], color="midnightblue",lwd=1.1)
psig =psig + theme_minimal() + xlab("snp") + ylab("beta")
cor_ref = data.frame(correlation = cor(Zscores[,which(grepl("Z",names(Zscores))), with=FALSE])[,1])
cor_ref["N"] = N_eff[row.names(cor_ref), "N_effective"]
print(head(cor_ref))
p_cor = ggplot(cor_ref, aes(x=N, y=correlation)) + geom_line() + geom_point(alpha=0.5) + theme_minimal()
beta_hat_prec = melt(Zscores[,c("V1","Beta_meta-analysis",paste0("Beta_",N_eff[N_eff$N_effective > 500,"Name"])), with=FALSE], id.vars=1)
psig_prec = ggplot(beta_hat_prec, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat_prec[variable=="Beta_meta-analysis",], color="midnightblue",lwd=1.1)
psig_prec = psig_prec + theme_minimal() + xlab("snp") + ylab("beta")
png("../../reports/figures/Signal_variability_simulated.png", width=900, height=325)
plot_grid(p_cor, psig, psig_prec, labels=c("A", "B", "C"), nrow=1)
dev.off()
beta_hat = melt(Zscores[,c(1,which(grepl("Beta_",names(Zscores)))), with=FALSE], id.vars=1)
psig = ggplot(beta_hat, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat[variable=="Beta_meta-analysis",], color="midnightblue",lwd=1.1)
psig =psig + theme_minimal() + xlab("snp") + ylab("beta")
panel_name = paste0("../../reports/figures/Signal_variability_simulated", cohort,"_",tag,".png")
cohort
if(cohort != '_small_cohort'){
beta_hat_prec = melt(Zscores[,c("V1","Beta_all_SNPs",paste0("Beta_",N_eff[N_eff$N_effective > 500,"study"])), with=FALSE], id.vars=1)
psig_prec = ggplot(beta_hat_prec, aes(x=V1, y=abs(value), group=variable)) + geom_line(alpha=0.2) + geom_line(data=beta_hat_prec[variable=="Beta_meta-analysis",], color="midnightblue",lwd=1.1)
psig_prec = psig_prec + theme_minimal() + xlab("snp") + ylab("beta")
panel_plot = plot_grid(p_cor, psig, psig_prec, labels=c("A", "B", "C"), nrow=1)
save_plot(panel_name, panel_plot,base_asp=1.8, bg='white' )
}else{
panel = plot_grid(p_cor, psig, labels=c("A", "B"), nrow=1)
panel
save_plot(panel_name, panel , bg='white' )
}
}
}
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment