Statistical plotting recipes for QTL analysis.
BigRiverQTLPlots.jl is a versatile plotting package built in the Julia programming language. The package consists of specific plotting recipes, designed to streamline data visualization and enhance the process of statistical analysis in genetic studies. It's particularly suited for QTL (Quantitative Trait Loci) and eQTL (expression Quantitative Trait Loci) analyses, providing key features for clear and intuitive data representation.
Here are the main functions provided by BigRiverQTLPlots.jl:
-
plot_QTL(): This function generates plots for LOD (logarithm of odds) scores with respect to marker positions. It's useful in viewing one genome scan result, or even multiple genome scan results, on a single plot. This helps to provide a broad overview of the QTL landscape. -
plot_eQTL(): This function is specifically designed to plot eQTL analysis results, assisting in visualization of genetic associations with gene expression levels.
To install BigRiverQTLPlots.jl, you can use Julia's package manager. Here is the command:
using Pkg
Pkg.add("BigRiverQTLPlots")After installing BigRiverQTLPlots.jl, you can include it in your Julia script using the following command:
using BigRiverQTLPlotsFrom there, you can start using plot_QTL and plot_eQTL to plot your data. For example:
# Assuming `single_results_perms` are restulting lod scores
# gInfo contains genotype information
plot_QTL(single_results_perms, gInfo)# Assuming `lod_scores` are in multipletraits_results
# pInfo contains phenotype information
# gInfo contains genotype information
# thresh is your LOD threshold value
plot_eQTL(multipletraits_results, pInfo, gInfo; threshold = 5.0)The following examples provide a basic idea of how to use the functions provided by BigRiverPlots.jl. Before proceeding, ensure that you've installed BigRiverPlots.jl and imported it into your Julia script using using BigRiverPlots along with Plots.jl. We use the package BulkLMM.jl to process the data.
using BigRiverQTLPlots, BulkLMM
using Plots
##############
# BXD spleen #
##############
########
# Data #
########
bulklmmdir = dirname(pathof(BulkLMM));
gmap_file = joinpath(bulklmmdir, "..", "data", "bxdData", "gmap.csv");
gInfo = BulkLMM.CSV.read(gmap_file, BulkLMM.DataFrames.DataFrame);
phenocovar_file = joinpath(bulklmmdir, "..", "data", "bxdData", "phenocovar.csv");
pInfo = BulkLMM.CSV.read(phenocovar_file, BulkLMM.DataFrames.DataFrame);
pheno_file = joinpath(bulklmmdir, "..", "data", "bxdData", "spleen-pheno-nomissing.csv");
pheno = BulkLMM.DelimitedFiles.readdlm(pheno_file, ',', header = false);
# exclude the header, the first (transcript ID)and the last columns (sex)
pheno_processed = pheno[2:end, 2:(end-1)] .* 1.0;
geno_file = joinpath(bulklmmdir, "..", "data", "bxdData", "spleen-bxd-genoprob.csv")
geno = BulkLMM.DelimitedFiles.readdlm(geno_file, ',', header = false);
geno_processed = geno[2:end, 1:2:end] .* 1.0;
#################
# Preprocessing #
#################
traitID = 1112;
pheno_y = reshape(pheno_processed[:, traitID], :, 1);
###########
# Kinship #
###########
kinship = calcKinship(geno_processed);
########
# Scan #
########
single_results_perms = scan(
pheno_y,
geno_processed,
kinship;
permutation_test = true,
nperms = 1000,
);
########
# Plot #
########
# to change bar color in the background use the keyword `barcolor`
plot_QTL(single_results_perms, gInfo, significance = [0.10, 0.05])
# or by specifying the thresholds values
thrs = BigRiverQTLPlots.perms_thresholds(single_results_perms.L_perms, [0.10, 0.05]);
plot_QTL(single_results_perms.lod, gInfo, thresholds = thrs);
The colors of the background bars and the threshold lines in the chart can be easily customized. To change the color of the background bars, use the barcolor argument, and to set the color of the threshold lines, use the thresholdcolor argument. These arguments accept standard Julia color codes, and HEX values.
Similarly, we can generate a manhattan plot using plot_manhattan function.
plot_manhattan(single_results_perms, gInfo, significance = [0.10, 0.05] );
Use the keywords manhattancolor to customize the scatter plot colors.
plot_manhattan(
single_results_perms, gInfo,
significance = [0.10, 0.05],
manhattancolor = ["#a6cee3", "#1f78b4"],
);
using BigRiverQTLPlots, BulkLMM
using Plots
##############
# BXD spleen #
##############
########
# Data #
########
bulklmmdir = dirname(pathof(BulkLMM));
gmap_file = joinpath(bulklmmdir, "..", "data", "bxdData", "gmap.csv");
gInfo = BulkLMM.CSV.read(gmap_file, BulkLMM.DataFrames.DataFrame);
phenocovar_file = joinpath(bulklmmdir, "..", "data", "bxdData", "phenocovar.csv");
pInfo = BulkLMM.CSV.read(phenocovar_file, BulkLMM.DataFrames.DataFrame);
pheno_file = joinpath(bulklmmdir, "..", "data", "bxdData", "spleen-pheno-nomissing.csv");
pheno = BulkLMM.DelimitedFiles.readdlm(pheno_file, ',', header = false);
geno_file = joinpath(bulklmmdir, "..", "data", "bxdData", "spleen-bxd-genoprob.csv")
geno = BulkLMM.DelimitedFiles.readdlm(geno_file, ',', header = false);
geno_processed = geno[2:end, 1:2:end] .* 1.0;
###########
# Kinship #
###########
kinship = calcKinship(geno_processed);
########
# Scan #
########
multipletraits_results, heritability_results = bulkscan_null(
pheno_processed,
geno_processed,
kinship,
)
########
# Plot #
########
plot_eQTL(multipletraits_results, pInfo, gInfo; threshold = 5.0);
Contributions to BigRiverQTLPlots.jl are welcome and appreciated. If you'd like to contribute, please fork the repository and make changes as you'd like. If you have any questions or issues, feel free to open an issue on the repository.
BigRiverQTLPlots.jl is licensed under the GNU AFFERO GENERAL PUBLIC LICENSE. For more information, please refer to the LICENSE file in the repository.
If you have any problems or questions using BigRiverQTLPlots.jl, please open an issue on the GitHub repository. We'll be happy to help!
React
Typescript
Django
Laravel
Facebook
Microsoft
Google
Alibaba
D3
Tencent