- Score computation workflow
- Getting started
GENOSCORES is a platform which encapsulates an R package, and a web platform which are linked by a public API. This provides flexibility of using GENOSCORES, as one is able to:
- Use R package with the example script directly (recommended).
- Write own R code to interact with the R package and/or the API.
- Write custom calls to the API in any language of choice.
The R package is by far the easiest way to use GENOSCORES to compute genotypic scores. One can simply change the parameters in the provided example script and run the analysis. The script will then call relevant package functions which will make API calls to retrieve the data automatically. Hence, no knowledge of the package, or the API is needed.
Score computation workflow
Genotypic scores are computed as a sum of genotypes, g1, weighted by the effect size estimate, β, and further multiplied by the inverse of the SNP-SNP correlation matrix2, R, to adjust the score for LD structure.
The score for an individual
tis given by a product of β for trait
tby the inverse of R by the vector of genotypes g for individual
Genotypes g are provided by the user in one of the supported formats: PLINK, dosage or bgendosage (bgen file transformed into dosage format with qctools). However, we only support input file preprocessing for files in PLINK format.
The weights (effect sizes) β are stored in
Effect sizes β for selected studies and/or trait types
are downloaded automatically during analysis. Summary statistics (weights) are
filtered at a p-value of
1E-5. This can be changed in
example script depending on the analysis needs.
We use a reference panel to compute the correlation matrix R. User is expected to provide a reference panel file in PLINK format. For example, we often use European samples in the 1000Genomes project as a reference panel. If the reference panel is not provided, LD adjustment will be omitted.
We compute genotypic scores separately for each trait-associated region of the genome. This means that for each GWAS requested from our database you may get more than 1 score computed.
Example: If the GWAS of gene expression for a given gene has detected a cis-eQTL and a trans-eQTL, then two scores would be computed corresponding to the cis locus and the trans locus.
Trait-associated regions are defined as genomic regions containing at least one
p-value < 1E-7. You can change this threshold in
example script if needed (e.g. use genome-wide significance).
The boundaries of each genomic region are defined by positions at which there
is a gap of at least one megabase from any other SNP in the filtered set of SNPs
for that GWAS. SNPs in the filtered set that are not assigned to trait-associated
regions are assigned to a global “background” score for the GWAS. The
"background" score has a score id with suffix
To compute a genome-wide score, you can simply add all the regional scores and the background score for a given GWAS.
You may need to run GENOSCORES on high performance computers depending on your use case4.
Create user account
Navigate to the registration page and create your user account.
This is needed to authenticate the requests to the web platform that R package makes on your behalf.
Download and install the package
Download the R package, and the example script.
To install the R package navigate to the containing directory and execute the following command:
R CMD INSTALL genoscores_xxx.tar.gz
xxx stands for the package version.
Install the dependencies
GENOSCORES R package requires the fallowing R packages to be installed:
"RMariaDB" "Rcpp" "bigmemory" "corrplot" "data.table" "doParallel" "foreach" "ggplot2" "gtools" "httr" "jsonlite" "Matrix" "RColorBrewer" "yaml" "BH" "RcppEigen"
To install the above package in a single command, simply run inside R session:
install.packages(c("RMariaDB", "Rcpp", "bigmemory", "corrplot", "data.table", "doParallel", "foreach", "ggplot2", "gtools", "httr", "jsonlite", "Matrix", "RColorBrewer", "yaml", "BH", "RcppEigen"))
Notes on installing bedtools
Pull the package with
Extract the package with
tar xzf bedtools-2.30.0.tar.gz
Compile and build bedtools:
Navigate to the source directory and run:
By default, bedtools will be installed in
/usr/local/bin. If this directory is not accessible, use a custom path:
make DESTDIR=/custom/dir install
do not forget to add that custom directory to your system path
Initialise the package
To use GENOSCORES you need to authorize the package to make web requests to the web application on your behalf. During initialization the package will also check that all the dependencies have been installed correctly.
Initialisation is interactive. In R session load GENOSCORES library and run initialisation like this:
This will start an interactive initialisation session. Follow the on-screen instructions.
Upon successful initialisation, a hidden file
.genoscores.yml will be
created in your
HOME directory. It contains credentials for authenticating
requests to the web platform. Please do not move,
copy or otherwise expose this information to any third parties.
At this point everything should be ready to run the analysis!
All the analysis input parameters are described in detail in the example script, and you are encouraged to examine this script carefully. Here we give a generic overview of all the configuration that the example script provides.
Specify which studies and/or trait types to use
You can select which GWAS to be used in the analysis by providing a vector of study ids and/or a vector of trait types. Check Studies page for details.
In example script this is controlled by the parameters:
trait.types or either of them
can be specified.
Configure system resources
This section of example script controls resource allocation and configures execution of the package functions in parallel mode.
Important parameters to tweak are:
Also consider parameter
memory.plink is the input file format is PLINK.
Set parameters for input data processing
The input file should contain genotype data from the patients for which the regional genotypic scores are to be computed. We currently support PLINK, dosage, bgendosage (bgen file transformed into dosage format with qctools) as input file formats. Note, that we currently support the input file pre-processing only for PLINK files.
For input data specified in PLINK format we offer two stages of data pre-processing:
- Update the SNP ids in the input file to the most recent build.
- Filter SNP ids in the input file to keep only those present in
the GENOSCORES database. Checkout
Updating of the genomic positions requires a liftOver chain file to lift positions to build hg38. If your data is already in build hg38 you can skip this step.
If you want to find out the build of your data, you can use
function shipped with the package.
Set parameters for score computation and adjustment for LD
This section of example script contains quite a few parameters which control various features of the analysis such as:
- splitting the genome into genomic regions
- filtering SNPs by
- LD adjustment reference panel
- LD adjustment methods
On linkage disequilibrium
The presence of linkage disequilibrium (LD) can affect the quality of the scores generated from the target input file, as the contribution of SNPs in high LD may overwhelm that of other SNPs. To avoid a potential distortion of the scores, the effect sizes (betas) of SNPs in high LD need to be adjusted before the computation of scores.
GENOSCORES is able to perform such LD adjustment automatically if a reference
population genome is provided in the
ld.refplinkfile option of
Therefore, it is strongly recommended to provide reference PLINK file for the
same population as in the target file.
As an example, below we present the list of steps required to generate such file using EUR ancestry in 1000 Genomes Project (other reference projects such as UK10K can be used as well):
- Download the reference genotypes from 1000 Genomes Project data portal. Make sure that the genotypes are in PLINK format and only samples of European ancestry are present, if your target population is of European ancestry.
- Update RSids to build hg38. You can use
prepare.plinkfile()from GENOSCORES to do this.
Now you can use this file as
ld.refplinkfile for LD adjustment in GENOSCORES.
When you execute the example script, GENOSCORES will compute a genotypic score for each sample in the target cohort (provided genotypes g) and each trait-associated region for the selected GWAS. See Score computation workflow for details about the score computation.
GENOSCORES will produce the following files and store them in the specified output directory.
Contains the (
number_of_samples X number_of_scores) matrix with the computed
The rows of the matrix contain IDs of the individuals. If the input file format is PLINK, then a composite ID5 is used.
The columns of the matrix contain score IDs. Consider score
1_0 as example.
The first digit
1 stands for
gwasid: identifier of a particular GWAS from
which β (weights) were used to compute the score.
The second digit
0 specifies that the score is a background score.
If the second digit is not
0, than the score is a regional score.
data.table containing aggregated information on scores, GWAS,
β (weights) used to compute these scores and studies
from which these weights come.
||ID of a GWAS used to compute a score.|
||ID of a genomic region of a score. If
||Chromosome that carries the top SNP of a score.|
||Genomic position of the top SNP of a score.|
||ID of a study containing the GWAS used in computing a score.|
||ID of a trait for which GWAS was performed.|
||ID of the top SNP of a score in DbSNP table.|
||Name of the trait for which GWAS was performed. For example
||Type of the trait. For example
||Tissue in which the trait was measured. For example
||The reference allele we have considered.|
||Effect size (weight) taken from the GWAS for a given trait which was used in computing a score.|
||Reference allele as reported by DbSNP|
||Alt allele as reported by DbSNP|
||PubMed identifier of a publication containing the GWAS used to compute a score.|
||Number of samples used in the study.|
||Number of cases if this is a case-control study.|
||Number of controls if this is a case-control study.|
||Genetic ancestry of the samples in the study.|
||The chromosome as an integer (with X=23, Y=24, etc).|
||SNP as provided by the user.|
||Genomic position of the beginning of the score.|
||Genomic position of the end of the score.|
||Minimum p-value recorded within a region.|
region makes a score ID.
data.table containing information for each score. Most columns here are
||Specifies how many SNPs contribute to the score.|
||Lists genes which are located close to the region of the score.|
Note that identical information is also saved to
For traits of type
QTLs, user may call function
annotate.qtls and pass
scoresinfo table loaded from
Annotate the table of scores information with gene positions of QTLs.
Currently, the trait types supported by this function are "transcript" and
"splicing". This procedure appends additional columns to
||Short name of a gene as specified in Ensembl genome browser|
||Is a distance from a gene to the score|
This file contains all the settings specified in the example script and used to compute scores (p-value thresholds, LD-adjustment parameters, etc). Full specification is available in example script. The key parameters are listed in Inputs. This is useful for reproducibility.
For example, if the input file format is PLINK, the genotypes are the scored (as 0, 1, 2)copies of the reference allele. ↩
We have implemented a pseudo-inverse solution to invert singular correlation matrices occurring in regions of high LD. ↩
We search for, download and process publicly available summary statistics from GWASs on a number of traits including: disease states (especially autoimmune diseases), clinical traits, gene expression in various tissues, metabolites and proteins, and others. You can check the Studies page for a complete list of studies and trait types currently available. ↩
GENOSCORES R package can be installed and run on LINUX and macOS. We do not currently support Windows. GENOSCORES performs expensive matrix operations on user computers to calculate scores. If large input datasets are used, it may be beneficial to use high performance multiple cores Linux machines or an HPC. ↩