| Title: | Colocalisation Tests of Two Genetic Traits |
|---|---|
| Description: | Performs the colocalisation tests described in Giambartolomei et al (2013) <doi:10.1371/journal.pgen.1004383>, Wallace (2020) <doi:10.1371/journal.pgen.1008720>, Wallace (2021) <doi:10.1371/journal.pgen.1009440>. |
| Authors: | Chris Wallace [aut, cre], Claudia Giambartolomei [aut], Vincent Plagnol [ctb] |
| Maintainer: | Chris Wallace <[email protected]> |
| License: | GPL |
| Version: | 5.2.3 |
| Built: | 2024-11-13 20:58:15 UTC |
| Source: | https://github.com/chr1swallace/coloc |
Performs the colocalisation tests described in Plagnol et al (2009) and Wallace et al (2020) and draws some plots.
Chris Wallace [email protected]
coloc functions need to be able to link summary stats from two different datasets and they do this through snp identifiers. This function takes the output of susie_rss() and adds snp identifiers. It is entirely the user's responsibility to ensure snp identifiers are in the correct order, coloc cannot make any sanity checks.
annotate_susie(res, snp, LD)annotate_susie(res, snp, LD)
res |
output of susie_rss() |
snp |
vector of snp identifiers |
LD |
matrix of LD (r) between snps in snp identifiers. Columns, rows should be named by a string that exists in the vector snp |
Note: this annotation step is not needed if you use runsusie() - this is only required if you use the susieR functions directly
res with column names added to some components
Chris Wallace
Internal function, approx.bf.estimates
approx.bf.estimates( z, V, type, suffix = NULL, sdY = 1, effect_priors = c(quant = 0.15, cc = 0.2) )approx.bf.estimates( z, V, type, suffix = NULL, sdY = 1, effect_priors = c(quant = 0.15, cc = 0.2) )
z |
normal deviate associated with regression coefficient and its variance |
V |
its variance |
type |
"quant" or "cc" |
suffix |
suffix to append to column names of returned data.frame |
sdY |
standard deviation of the trait. If not supplied, will be estimated. |
Calculate approximate Bayes Factors using supplied variance of the regression coefficients
data.frame containing lABF and intermediate calculations
Vincent Plagnol, Chris Wallace
Internal function, approx.bf.p
approx.bf.p(p, f, type, N, s, suffix = NULL)approx.bf.p(p, f, type, N, s, suffix = NULL)
p |
p value |
f |
MAF |
type |
"quant" or "cc" |
N |
sample size |
s |
proportion of samples that are cases, ignored if type=="quant" |
suffix |
suffix to append to column names of returned data.frame |
Calculate approximate Bayes Factors
data.frame containing lABF and intermediate calculations
Claudia Giambartolomei, Chris Wallace
Convert binomial to linear regression
bin2lin(D, doplot = FALSE)bin2lin(D, doplot = FALSE)
D |
standard format coloc dataset |
doplot |
plot results if TRUE - useful for debugging |
Estimate beta and varbeta if a linear regression had been run on a binary outcome, given log OR and their variance + MAF in controls
sets beta = cov(x,y)/var(x) varbeta = (var(y)/var(x) - cov(x,y)^2/var(x)^2)/N
D, with original beta and varbeta in beta.bin, varbeta.bin, and beta and varbeta updated to linear estimates
Chris Wallace
check alignment between beta and LD
check_alignment(D, thr = 0.2, do_plot = TRUE) check.alignment(...)check_alignment(D, thr = 0.2, do_plot = TRUE) check.alignment(...)
D |
a coloc dataset |
thr |
plot SNP pairs in absolute LD > thr |
do_plot |
if TRUE (default) plot the diagnostic |
... |
arguments passed to check_alignment() |
proportion of pairs that are positive
Chris Wallace
Check coloc dataset inputs for errors
check_dataset(d, suffix = "", req = c("type", "snp"), warn.minp = 1e-06) check.dataset(...)check_dataset(d, suffix = "", req = c("type", "snp"), warn.minp = 1e-06) check.dataset(...)
d |
dataset to check |
suffix |
string to identify which dataset (1 or 2) |
req |
names of elements that must be present |
warn.minp |
print warning if no p value < warn.minp |
... |
arguments passed to check_dataset() |
A coloc dataset is a list, containing a mixture of vectors capturing quantities that vary between snps (these vectors must all have equal length) and scalars capturing quantities that describe the dataset.
Coloc is flexible, requiring perhaps only p values, or z scores, or effect estimates and standard errors, but with this flexibility, also comes difficulties describing exactly the combinations of items required.
Required vectors are some subset of
regression coefficient for each SNP from dataset 1
variance of beta
P-values for each SNP in dataset 1
minor allele frequency of the variants
a character vector of snp ids, optional. It will be used to merge dataset1 and dataset2 and will be retained in the results.
Preferably, give beta and varbeta. But if these are not available, sufficient statistics can be approximated from pvalues and MAF.
Required scalars are some subset of
Number of samples in dataset 1
the type of data in dataset 1 - either "quant" or "cc" to denote quantitative or case-control
for a case control dataset, the proportion of samples in dataset 1 that are cases
for a quantitative trait, the population standard deviation of the trait. if not given, it can be estimated from the vectors of varbeta and MAF
You must always give type. Then,
type=="cc"s
type=="quant" and sdY knownsdY
N
If sdY is unknown, it will be approximated, and this will require
sdY
beta, varbeta, N, MAF
Optional vectors are
a vector of snp positions, required for plot_dataset
check_dataset calls stop() unless a series of expectations on dataset
input format are met
This is a helper function for use by other coloc functions, but you can use it directly to check the format of a dataset to be supplied to coloc.abf(), coloc.signals(), finemap.abf(), or finemap.signals().
NULL if no errors found
Chris Wallace
Simulated data to use in testing and vignettes in the coloc package
data(coloc_test_data)data(coloc_test_data)
A four of two coloc-style datasets. Elements D1 and D2 have a single shared causal variant, and 50 SNPs. Elements D3 and D4 have 100 SNPs, one shared causal variant, and one variant unique to D3. Use these as examples of what a coloc-style dataset for a quantitative trait should look like.
data(coloc_test_data) names(coloc_test_data) str(coloc_test_data$D1) check_dataset(coloc_test_data$D1) # should return NULL if data structure is okdata(coloc_test_data) names(coloc_test_data) str(coloc_test_data$D1) check_dataset(coloc_test_data$D1) # should return NULL if data structure is ok
Bayesian colocalisation analysis
coloc.abf( dataset1, dataset2, MAF = NULL, p1 = 1e-04, p2 = 1e-04, p12 = 1e-05, ... )coloc.abf( dataset1, dataset2, MAF = NULL, p1 = 1e-04, p2 = 1e-04, p12 = 1e-05, ... )
dataset1 |
a list with specifically named elements defining
the dataset to be analysed. See |
dataset2 |
as above, for dataset 2 |
MAF |
Common minor allele frequency vector to be used for both dataset1 and dataset2, a shorthand for supplying the same vector as parts of both datasets |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
... |
used to pass parameters to approx.bf.estimates, in particular the effect_priors parameter |
This function calculates posterior probabilities of different causal variant configurations under the assumption of a single causal variant for each trait.
If regression coefficients and variances are available, it calculates Bayes factors for association at each SNP. If only p values are available, it uses an approximation that depends on the SNP's MAF and ignores any uncertainty in imputation. Regression coefficients should be used if available.
a list of two data.frames:
summary is a vector giving the number of SNPs analysed, and the posterior probabilities of H0 (no causal variant), H1 (causal variant for trait 1 only), H2 (causal variant for trait 2 only), H3 (two distinct causal variants) and H4 (one common causal variant)
results is an annotated version of the input data containing log Approximate Bayes Factors and intermediate calculations, and the posterior probability SNP.PP.H4 of the SNP being causal for the shared signal if H4 is true. This is only relevant if the posterior support for H4 in summary is convincing.
Claudia Giambartolomei, Chris Wallace
Colocalise two datasets represented by Bayes factors
coloc.bf_bf( bf1, bf2, p1 = 1e-04, p2 = 1e-04, p12 = 5e-06, overlap.min = 0.5, trim_by_posterior = TRUE )coloc.bf_bf( bf1, bf2, p1 = 1e-04, p2 = 1e-04, p12 = 5e-06, overlap.min = 0.5, trim_by_posterior = TRUE )
bf1 |
named vector of log BF, or matrix of BF with colnames (cols=snps, rows=signals) |
bf2 |
named vector of log BF, or matrix of BF with colnames (cols=snps, rows=signals) |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
overlap.min |
see trim_by_posterior |
trim_by_posterior |
it is important that the signals to be colocalised are covered by adequate numbers of snps in both datasets. If TRUE, signals for which snps in common do not capture least overlap.min proportion of their posteriors support are dropped and colocalisation not attempted. |
This is the workhorse behind many coloc functions
coloc.signals style result
Chris Wallace
Bayesian colocalisation analysis, detailed output
coloc.detail( dataset1, dataset2, MAF = NULL, p1 = 1e-04, p2 = 1e-04, p12 = 1e-05 )coloc.detail( dataset1, dataset2, MAF = NULL, p1 = 1e-04, p2 = 1e-04, p12 = 1e-05 )
dataset1 |
a list with specifically named elements defining
the dataset to be analysed. See |
dataset2 |
as above, for dataset 2 |
MAF |
Common minor allele frequency vector to be used for both dataset1 and dataset2, a shorthand for supplying the same vector as parts of both datasets |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
This function replicates coloc.abf, but outputs more detail for further processing using coloc.process
Intended to be called internally by coloc.signals
a list of three data.tabless:
summary is a vector giving the number of SNPs analysed, and the posterior probabilities of H0 (no causal variant), H1 (causal variant for trait 1 only), H2 (causal variant for trait 2 only), H3 (two distinct causal variants) and H4 (one common causal variant)
df is an annotated version of the input data containing log Approximate Bayes Factors and intermediate calculations, and the posterior probability SNP.PP.H4 of the SNP being causal for the shared signal
df3 is the same for all 2 SNP H3 models
Chris Wallace
Internal helper function
coloc.process( obj, hits1 = NULL, hits2 = NULL, LD = NULL, r2thr = 0.01, p1 = 1e-04, p2 = 1e-04, p12 = 1e-06, LD1 = LD, LD2 = LD, mode = c("iterative", "allbutone") )coloc.process( obj, hits1 = NULL, hits2 = NULL, LD = NULL, r2thr = 0.01, p1 = 1e-04, p2 = 1e-04, p12 = 1e-06, LD1 = LD, LD2 = LD, mode = c("iterative", "allbutone") )
obj |
object returned by coloc.detail() |
hits1 |
lead snps for trait 1. If length > 1, will use masking |
hits2 |
lead snps for trait 2. If length > 1, will use masking |
LD |
named LD matrix (for masking) |
r2thr |
r2 threshold at which to mask |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
LD1 |
named LD matrix (for masking) for trait 1 only |
LD2 |
named LD matrix (for masking) for trait 2 only |
mode |
either "iterative" (default) - successively condition on signals or "allbutone" - find all putative signals and condition on all but one of them in each analysis |
data.table of coloc results
Chris Wallace
New coloc function, builds on coloc.abf() by allowing for multiple independent causal variants per trait through conditioning or masking.
coloc.signals( dataset1, dataset2, MAF = NULL, LD = NULL, method = c("single", "cond", "mask"), mode = c("iterative", "allbutone"), p1 = 1e-04, p2 = 1e-04, p12 = NULL, maxhits = 3, r2thr = 0.01, pthr = 1e-06 )coloc.signals( dataset1, dataset2, MAF = NULL, LD = NULL, method = c("single", "cond", "mask"), mode = c("iterative", "allbutone"), p1 = 1e-04, p2 = 1e-04, p12 = NULL, maxhits = 3, r2thr = 0.01, pthr = 1e-06 )
dataset1 |
a list with specifically named elements defining
the dataset to be analysed. See |
dataset2 |
as above, for dataset 2 |
MAF |
Common minor allele frequency vector to be used for both dataset1 and dataset2, a shorthand for supplying the same vector as parts of both datasets |
LD |
required if method="cond". matrix of genotype correlation (ie r, not r^2) between SNPs. If dataset1 and dataset2 may have different LD, you can instead add LD=LD1 to the list of dataset1 and a different LD matrix for dataset2 |
method |
default "" means do no conditioning, should return similar to coloc.abf. if method="cond", then use conditioning to coloc multiple signals. if method="mask", use masking to coloc multiple signals. if different datasets need different methods (eg LD is only available for one of them) you can set method on a per-dataset basis by adding method="..." to the list for that dataset. |
mode |
"iterative" or "allbutone". Easiest understood with an example. Suppose there are 3 signal SNPs detected for trait 1, A, B, C and only one for trait 2, D. Under "iterative" mode, 3 coloc will be performed: * trait 1 - trait 2 * trait 1 conditioned on A - trait 2 * trait 1 conditioned on A+B - trait 2 Under "allbutone" mode, they would be * trait 1 conditioned on B+C - trait 2 * trait 1 conditioned on A+C - trait 2 * trait 1 conditioned on A+B - trait 2 Only iterative mode is supported for method="mask". The allbutone mode is optimal if the signals are known with certainty (which they never are), because it allows each signal to be tested without influence of the others. When there is uncertainty, it may make sense to use iterative mode, because the strongest signals aren't affected by conditioning incorrectly on weaker secondary and less certain signals. |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
maxhits |
maximum number of levels to condition/mask |
r2thr |
if masking, the threshold on r2 should be used to call two signals independent. our experience is that this needs to be set low to avoid double calling the same strong signal. |
pthr |
if masking or conditioning, what p value threshold to call a secondary hit "significant" |
data.table of coloc results, one row per pair of lead snps detected in each dataset
Chris Wallace
colocalisation with multiple causal variants via SuSiE
coloc.susie( dataset1, dataset2, back_calculate_lbf = FALSE, susie.args = list(), ... )coloc.susie( dataset1, dataset2, back_calculate_lbf = FALSE, susie.args = list(), ... )
dataset1 |
either a coloc-style input dataset (see check_dataset), or the result of running runsusie on such a dataset |
dataset2 |
either a coloc-style input dataset (see check_dataset), or the result of running runsusie on such a dataset |
back_calculate_lbf |
by default, use the log Bayes factors returned by susie_rss. It is also possible to back-calculate these from the posterior probabilities. It is not advised to set this to TRUE, the option exists really for testing purposes only. |
susie.args |
a named list of additional arguments to be passed to runsusie |
... |
other arguments passed to coloc.bf_bf, in particular prior values for causal association with one trait (p1, p2) or both (p12) |
a list, containing elements * summary a data.table of posterior probabilities of each global hypothesis, one row per pairwise comparison of signals from the two traits * results a data.table of detailed results giving the posterior probability for each snp to be jointly causal for both traits assuming H4 is true. Please ignore this column if the corresponding posterior support for H4 is not high. * priors a vector of the priors used for the analysis
Chris Wallace
coloc for susie output + a separate BF matrix
coloc.susie_bf( dataset1, bf2, p1 = 1e-04, p2 = 1e-04, p12 = 5e-06, susie.args = list(), ... )coloc.susie_bf( dataset1, bf2, p1 = 1e-04, p2 = 1e-04, p12 = 5e-06, susie.args = list(), ... )
dataset1 |
a list with specifically named elements defining
the dataset to be analysed. See |
bf2 |
named vector of log BF, names are snp ids and will be matched to column names of susie object's alpha |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
susie.args |
named list of arguments to be passed to susieR::susie_rss() |
... |
other arguments passed to coloc.bf_bf, in particular prior values for causal association with one trait (p1, p2) or both (p12) |
coloc.signals style result
Chris Wallace
Internal function, calculate posterior probabilities for configurations, given logABFs for each SNP and prior probs
combine.abf(l1, l2, p1, p2, p12, quiet = FALSE)combine.abf(l1, l2, p1, p2, p12, quiet = FALSE)
l1 |
merged.df$lABF.df1 |
l2 |
merged.df$lABF.df2 |
p1 |
prior probability a SNP is associated with trait 1, default 1e-4 |
p2 |
prior probability a SNP is associated with trait 2, default 1e-4 |
p12 |
prior probability a SNP is associated with both traits, default 1e-5 |
quiet |
don't print posterior summary if TRUE. default=FALSE |
named numeric vector of posterior probabilities
Claudia Giambartolomei, Chris Wallace
Internal helper function for est_all_cond
est_cond(x, LD, YY, sigsnps, xtx = NULL)est_cond(x, LD, YY, sigsnps, xtx = NULL)
x |
coloc dataset |
LD |
named matrix of r |
YY |
sum((Y-Ybar)^2) |
sigsnps |
names of snps to jointly condition on |
xtx |
optional, matrix X'X where X is the genotype matrix. If not available, will be estimated from LD, MAF, beta and sample size (the last three should be part of the coloc dataset) |
data.table giving snp, beta and varbeta on remaining snps after conditioning
Chris Wallace
Estimate single snp frequency distibutions
estgeno.1.ctl(f) estgeno.1.cse(G0, b)estgeno.1.ctl(f) estgeno.1.cse(G0, b)
f |
MAF |
G0 |
single snp frequency in controls (vector of length 3) - obtained from estgeno.1.ctl |
b |
log odds ratio |
relative frequency of genotypes 0, 1, 2
Chris Wallace
estgeno2
Internal helper function
find.best.signal(D)find.best.signal(D)
D |
standard format coloc dataset |
z at most significant snp, named by that snp id
Chris Wallace
tries to be smart about detecting the interesting subregion to finemap/coloc.
findends(d, maxz = 4, maxr2 = 0.1, do.plot = FALSE)findends(d, maxz = 4, maxr2 = 0.1, do.plot = FALSE)
d |
a coloc dataset |
maxz |
keep all snps between the leftmost and rightmost snp with |z| > maxz |
maxr2 |
expand window to keep all snps between snps with r2 > maxr2 with the left/rightmost snps defined by the maxz threshold |
do.plot |
if TRUE, plot dataset + boundaries |
logical vector of length d$position indicating which snps to keep
Chris Wallace
findpeaks
tries to be smart about detecting the interesting subregion to finemap/coloc.
findpeaks(d, maxz = 4, maxr2 = 0.1, do.plot = FALSE)findpeaks(d, maxz = 4, maxr2 = 0.1, do.plot = FALSE)
d |
a coloc dataset |
maxz |
keep all snps between the leftmost and rightmost snp with |z| > maxz |
maxr2 |
expand window to keep all snps between snps with r2 > maxr2 with the left/rightmost snps defined by the maxz threshold |
do.plot |
if TRUE, plot dataset + boundaries |
Differs from findends by finding multiple separate regions if there are multiple peaks
logical vector of length d$position indicating which snps to keep
Chris Wallace
findends
Bayesian finemapping analysis
finemap.abf(dataset, p1 = 1e-04)finemap.abf(dataset, p1 = 1e-04)
dataset |
a list with specifically named elements defining the dataset
to be analysed. See |
p1 |
prior probability a SNP is associated with the trait 1, default 1e-4 |
This function calculates posterior probabilities of different causal variant for a single trait.
If regression coefficients and variances are available, it calculates Bayes factors for association at each SNP. If only p values are available, it uses an approximation that depends on the SNP's MAF and ignores any uncertainty in imputation. Regression coefficients should be used if available.
a data.frame:
an annotated version of the input data containing log Approximate Bayes Factors and intermediate calculations, and the posterior probability of the SNP being causal
Chris Wallace
Finemap one dataset represented by Bayes factors
finemap.bf(bf1, p1 = 1e-04)finemap.bf(bf1, p1 = 1e-04)
bf1 |
named vector of log BF, or matrix of log BF with colnames (cols=snps, rows=signals) |
p1 |
prior probability a SNP is associated with the trait 1, default 1e-4 |
This is the workhorse behind many finemap functions
finemap.signals style result
Chris Wallace
This is an analogue to finemap.abf, adapted to find multiple signals where they exist, via conditioning or masking - ie a stepwise procedure
finemap.signals( D, LD = D$LD, method = c("single", "mask", "cond"), r2thr = 0.01, sigsnps = NULL, pthr = 1e-06, maxhits = 3, return.pp = FALSE )finemap.signals( D, LD = D$LD, method = c("single", "mask", "cond"), r2thr = 0.01, sigsnps = NULL, pthr = 1e-06, maxhits = 3, return.pp = FALSE )
D |
list of summary stats for a single disease, see check_dataset |
LD |
matrix of signed r values (not rsq!) giving correlation between SNPs |
method |
if method="cond", then use conditioning to coloc multiple signals. The default is mask - this is less powerful, but safer because it does not assume that the LD matrix is properly allelically aligned to estimated effect |
r2thr |
if mask==TRUE, all snps will be masked with r2 > r2thr with any sigsnps. Otherwise ignored |
sigsnps |
SNPs already deemed significant, to condition on or mask, expressed as a numeric vector, whose names are the snp names |
pthr |
when p > pthr, stop successive searching |
maxhits |
maximum depth of conditioning. procedure will stop if p > pthr OR abs(z)<zthr OR maxhits hits have been found. |
return.pp |
if FALSE (default), just return the hits. Otherwise return vectors of PP |
mask |
use masking if TRUE, otherwise conditioning. defaults to TRUE |
list of successively significant fine mapped SNPs, named by the SNPs
Chris Wallace
generic convenience function to convert logbf matrix to PP matrix
logbf_to_pp(bf, pi, last_is_null)logbf_to_pp(bf, pi, last_is_null)
bf |
an L by p or p+1 matrix of log Bayes factors |
pi |
either a scalar representing the prior probability for any snp to be causal, or a full vector of per snp / null prior probabilities |
last_is_null |
TRUE if last value of the bf vector or last column of a bf matrix relates to the null hypothesis of no association. This is standard for SuSiE results, but may not be for BF constructed in other ways. |
matrix of posterior probabilities, same dimensions as bf
Chris Wallace
Internal function, logdiff
logdiff(x, y)logdiff(x, y)
x |
numeric |
y |
numeric |
This function calculates the log of the difference of the exponentiated logs taking out the max, i.e. insuring that the difference is not negative
max(x) + log(exp(x - max(x,y)) - exp(y-max(x,y)))
Chris Wallace
Internal function, logsum
logsum(x)logsum(x)
x |
numeric vector |
This function calculates the log of the sum of the exponentiated logs taking out the max, i.e. insuring that the sum is not Inf
max(x) + log(sum(exp(x - max(x))))
Claudia Giambartolomei
Internal helper function for finemap.signals
map_cond(D, LD, YY, sigsnps = NULL)map_cond(D, LD, YY, sigsnps = NULL)
D |
dataset in standard coloc format |
LD |
named matrix of r |
YY |
sum(y^2) |
sigsnps |
names of snps to mask |
named numeric - Z score named by snp
Chris Wallace
Internal helper function for finemap.signals
map_mask(D, LD, r2thr = 0.01, sigsnps = NULL)map_mask(D, LD, r2thr = 0.01, sigsnps = NULL)
D |
dataset in standard coloc format |
LD |
named matrix of r |
r2thr |
mask all snps with r2 > r2thr with any in sigsnps |
sigsnps |
names of snps to mask |
named numeric - Z score named by snp
Chris Wallace
Plot a coloc structured dataset
plot_dataset( d, susie_obj = NULL, highlight_list = NULL, alty = NULL, ylab = "-log10(p)", show_legend = TRUE, color = c("dodgerblue2", "green4", "#6A3D9A", "#FF7F00", "gold1", "skyblue2", "#FB9A99", "palegreen2", "#CAB2D6", "#FDBF6F", "gray70", "khaki2", "maroon", "orchid1", "deeppink1", "blue1", "steelblue4", "darkturquoise", "green1", "yellow4", "yellow3", "darkorange4", "brown"), ... ) plot_dataset( d, susie_obj = NULL, highlight_list = NULL, alty = NULL, ylab = "-log10(p)", show_legend = TRUE, color = c("dodgerblue2", "green4", "#6A3D9A", "#FF7F00", "gold1", "skyblue2", "#FB9A99", "palegreen2", "#CAB2D6", "#FDBF6F", "gray70", "khaki2", "maroon", "orchid1", "deeppink1", "blue1", "steelblue4", "darkturquoise", "green1", "yellow4", "yellow3", "darkorange4", "brown"), ... )plot_dataset( d, susie_obj = NULL, highlight_list = NULL, alty = NULL, ylab = "-log10(p)", show_legend = TRUE, color = c("dodgerblue2", "green4", "#6A3D9A", "#FF7F00", "gold1", "skyblue2", "#FB9A99", "palegreen2", "#CAB2D6", "#FDBF6F", "gray70", "khaki2", "maroon", "orchid1", "deeppink1", "blue1", "steelblue4", "darkturquoise", "green1", "yellow4", "yellow3", "darkorange4", "brown"), ... ) plot_dataset( d, susie_obj = NULL, highlight_list = NULL, alty = NULL, ylab = "-log10(p)", show_legend = TRUE, color = c("dodgerblue2", "green4", "#6A3D9A", "#FF7F00", "gold1", "skyblue2", "#FB9A99", "palegreen2", "#CAB2D6", "#FDBF6F", "gray70", "khaki2", "maroon", "orchid1", "deeppink1", "blue1", "steelblue4", "darkturquoise", "green1", "yellow4", "yellow3", "darkorange4", "brown"), ... )
d |
a coloc dataset |
susie_obj |
optional, the output of a call to runsusie() |
highlight_list |
optional, a list of character vectors. any snp in the character vector will be highlighted, using a different colour for each list. |
alty |
default is to plot a standard manhattan. If you wish to plot a different y value, pass it here. You may also want to change ylab to describe what you are plotting. |
ylab |
label for y axis, default is -log10(p) and assumes you are plotting a manhattan |
show_legend |
optional, show the legend or not. default is TRUE |
color |
optional, specify the colours to use for each credible set when susie_obj is supplied. Default is shamelessly copied from susieR::susie_plot() so that colours will match |
... |
other arguments passed to the base graphics plot() function |
Chris Wallace
plot a coloc_abf object
## S3 method for class 'coloc_abf' plot(x, ...)## S3 method for class 'coloc_abf' plot(x, ...)
x |
coloc_abf object to be plotted |
... |
other arguments |
ggplot object
Chris Wallace
Print summary of a coloc.abf run
## S3 method for class 'coloc_abf' print(x, ...)## S3 method for class 'coloc_abf' print(x, ...)
x |
object of class |
... |
optional arguments: "trait1" name of trait 1, "trait2" name of trait 2 |
x, invisibly
Chris Wallace
Internal function, process each dataset list for coloc.abf.
process.dataset(d, suffix, ...)process.dataset(d, suffix, ...)
d |
list |
suffix |
"df1" or "df2" |
... |
used to pass parameters to approx.bf.estimates, in particular the effect_priors parameter |
Made public for another package to use, but not intended for users to use.
data.frame with log(abf) or log(bf)
Chris Wallace
run susie_rss storing some additional information for coloc
runsusie( d, suffix = 1, maxit = 100, repeat_until_convergence = TRUE, s_init = NULL, ... )runsusie( d, suffix = 1, maxit = 100, repeat_until_convergence = TRUE, s_init = NULL, ... )
d |
coloc dataset, must include LD (signed correlation matrix) and N (sample size) |
suffix |
suffix label that will be printed with any error messages |
maxit |
maximum number of iterations for the first run of susie_rss(). If susie_rss() does not report convergence, runs will be extended assuming repeat_until_convergence=TRUE. Most users will not need to change this default. |
repeat_until_convergence |
keep running until susie_rss() indicates convergence. Default TRUE. If FALSE, susie_rss() will run with maxit iterations, and if not converged, runsusie() will error. Most users will not need to change this default. |
s_init |
used internally to extend runs that haven't converged. don't use. |
... |
arguments passed to susie_rss. In particular, if you want to match some coloc defaults, set
otherwise susie_rss will estimate the prior variance itself |
results of a susie_rss run, with some added dimnames
Chris Wallace
library(coloc) data(coloc_test_data) result=runsusie(coloc_test_data$D1) summary(result)library(coloc) data(coloc_test_data) result=runsusie(coloc_test_data$D1) summary(result)
Estimate trait standard deviation given vectors of variance of coefficients, MAF and sample size
sdY.est(vbeta, maf, n)sdY.est(vbeta, maf, n)
vbeta |
vector of variance of coefficients |
maf |
vector of MAF (same length as vbeta) |
n |
sample size |
Estimate is based on var(beta-hat) = var(Y) / (n * var(X)) var(X) = 2maf(1-maf) so we can estimate var(Y) by regressing n*var(X) against 1/var(beta)
estimated standard deviation of Y
Chris Wallace
Shows how prior and posterior per-hypothesis probabilities change as a function of p12
sensitivity( obj, rule = "", dataset1 = NULL, dataset2 = NULL, npoints = 100, doplot = TRUE, plot.manhattans = TRUE, preserve.par = FALSE, row = 1 )sensitivity( obj, rule = "", dataset1 = NULL, dataset2 = NULL, npoints = 100, doplot = TRUE, plot.manhattans = TRUE, preserve.par = FALSE, row = 1 )
obj |
output of coloc.detail or coloc.process |
rule |
a decision rule. This states what values of posterior probabilities "pass" some threshold. This is a string which will be parsed and evaluated, better explained by examples. "H4 > 0.5" says post prob of H4 > 0.5 is a pass. "H4 > 0.9 & H4/H3 > 3" says post prob of H4 must be > 0.9 AND it must be at least 3 times the post prob of H3." |
dataset1 |
optional the dataset1 used to run SuSiE. This will be used to make a Manhattan plot if plot.manhattans=TRUE. |
dataset2 |
optional the dataset2 used to run SuSiE. This will be used to make a Manhattan plot if plot.manhattans=TRUE. |
npoints |
the number of points over which to evaluate the prior values for p12, equally spaced on a log scale between p1*p2 and min(p1,p2) - these are logical limits on p12, but not scientifically sensible values. |
doplot |
draw the plot. set to FALSE if you want to just evaluate the prior and posterior matrices and work with them yourself |
plot.manhattans |
if TRUE, show Manhattans of input data |
preserve.par |
if TRUE, do not change par() of current graphics device - this is to allow sensitivity plots to be incoporated into a larger set of plots, or to be plot one per page on a pdf, for example |
row |
when coloc.signals() has been used and multiple rows are returned in the coloc summary, which row to plot |
Function is called mainly for plotting side effect. It draws two plots, showing how prior and posterior probabilities of each coloc hypothesis change with changing p12. A decision rule sets the values of the posterior probabilities considered acceptable, and is used to shade in green the region of the plot for which the p12 prior would give and acceptable result. The user is encouraged to consider carefully whether some prior values shown within the green shaded region are sensible before accepting the hypothesis. If no shading is shown, then no priors give rise to an accepted result.
list of 3: prior matrix, posterior matrix, and a pass/fail indicator (returned invisibly)
Chris Wallace
Subset a coloc dataset
subset_dataset(dataset, index)subset_dataset(dataset, index)
dataset |
coloc dataset |
index |
vector of indices of snps to KEEP |
a copy of dataset, with only the data relating to snps in index remaining
Chris Wallace
variance of MLE of beta for quantitative trait, assuming var(y)=1
Var.data(f, N)Var.data(f, N)
f |
minor allele freq |
N |
sample number |
Internal function
variance of MLE beta
Claudia Giambartolomei
variance of MLE of beta for case-control
Var.data.cc(f, N, s)Var.data.cc(f, N, s)
f |
minor allele freq |
N |
sample number |
s |
??? |
Internal function
variance of MLE beta
Claudia Giambartolomei