evian_logit {evian}R Documentation

Evidential Analysis for dichotomous outcome data using logistic regression models

Description

Calculate and plot the likelihood intervals for genetic association in a genomic region of interest with a dichotomous trait. Covariates can be accommodated.

Usage

evian_logit(data, map, ycol, xcols, formula, robust = FALSE,
            model = 'all', rfdum = 'alphanum', m = 200, bse, 
            lolim = log(0.025), hilim = log(4), 
            plot_li = TRUE, yaxislim = 10, showmaxlr = 3, 
            kcutoff = 32, grays = TRUE, save_as_png = FALSE)

Arguments

data a data frame whose first six columns are in standard linkage format: FID, IID, PAT, MAT, SEX (1 = M, 2 = F), PHENOTYPE (0,1), followed by columns of non-genetic covariates, and then the columns of genotype data coded (0,1,2,NA) as produced in the ped file by the plink software with option --recodeA.
map a data frame with three column header: chr, snp, pos, ie. map file in plink format.
ycol numeric, the column number in the data data frame for the Y outcome variable, ie. PHENOTYPE or caseflag.
xcols range, the column range in the data data frame for the snp data, ie. 8:17.
formula string containing the glm formula for Y, e.g., Y ~ x + Z\$age for models where the snp covariate is a 1df parameterization, or Y ~ x + x1 + Z\$age for 2df snp parameterization. Use as.factor(Z\$var) if the variable is categorical. Note, the formula must be at least Y ~ x for 1df, and Y ~ x + x1 for 2df.
robust logical, if TRUE, then a robust adjustment is applied to the likelihood function to account for the cluster nature in the data, e.g. family id, FID.
model a string specifying the mode of inheritance parameterization:
additive, dominant, recessive, overdominance, 2df or all.
rfdum string to designate the reference group for categorical variables. Options: alphanum orders levels alphanumerically, most chooses the level observed most frequently, least chooses the level observed least frequently.
m numeric, the density of the grid at which to compute the standardized likelihood function. A beta grid is defined as the grid of values for the snp parameter at which to evaluate the likelihood function.
bse numeric, the number of beta standard errors to utilize in constraining the beta grid limits, e.g., beta grid is evaluated at hat{β} +/- 5 s.e.
lolim numeric, instead of the number of s.e., specify the lower limit for the grid, or the minimum value of the log(OR), at which to calculate the likelihood function.
hilim numeric, instead of the number of s.e., specify the upper limit for the grid, or the maximum value of the log(OR), at which to calculate the likelihood funciton.
plot_li logical, produce the likelihood intervals plot in graphics window (X11 or Quartz).
yaxislim y axis limit for the odds ratio, default = 10.
showmaxlr specifies the number of SNPs for which the user wants the maximum likelihood ratios displayed in text on the plot, default = 3.
kcutoff the strength of evidence criterion k: 32, 100, 1000, default = 32.
grays logical, if TRUE, then individual intervals for a given SNP will be in gray if the 1/k likelihood interval includes OR = 1 as a plausible value.
save_as_png logical, if TRUE the likelihood intervals plot will be saved as a png file. Filename will be a time stamp in current directory.

Details

logit_snp is main function called to calculate the 1/k likelihood intervals for the additive, dominant, recessive, overdominance and 2df genotypic models.

Usually the data inputs, ie. data and map, are conformed to the plink format when run with the --recodeA option.

The formula must be at least Y ~ x for snp 1df models, and Y ~ x + x1 for snp 2df models. Exact variable syntax and character spacing in formula is important for correct parsing and computation. In other words, a single whitespace separates each variable or symbol in formula.

Value

This function produces a data frame named according to the selected model, ie. data.additive, with the following 11 column headings:
snp, pos, maxlr, mle, k32_low, k32_hi, k100_low, k100_hi, k1000_low, k1000_hi, flip
When the value of mle is less than one, the k intervals are inverted and flip=1, otherwise flip=0.
By default, the likelihood intervals are plotted in an X11 window, and/or the image is saved in png file if save_as_png is TRUE. In the plot, the SNPs with OR = 1 in the 1/k likelihood interval will be in gray if grays = TRUE. Those with OR >= 1 not in the 1/k likelihood interval will be in color with the 1/32 intervals in green, 1/100 in red and 1/1000 in blue. maxLR is defined as

L(hat{theta})/L(1)

, where hat{theta} is the maximum likelihood estimate of the log(OR).
Rename the data frame in your R workspace before running another analysis to prevent overwriting results.

Note

When the number of beta standard errors bse is NOT defined, then the boundaries of the beta grid will use the default lolim/hilim values, or those defined by the user. However, if bse is defined, then lolim/hilim values will be ignored.

In some cases the beta grid (using bse or lolim/hilim), may need to be increased substantially (bse as large as 15) if covariates are present in the formula.

Depending on the size of the dataset, this function may require long compute times. A typical runtime for 1df models, ie. additive, dominant or recessive, without any covariates on a single processor of a 2.53 GHz Intel Core 2 Duo will be around 3 minutes for a dataset containing 50 snps with 250 individuals.

Finally, estimation may become inaccurate with large number of correlated covariates, similar to known limitation of profile likelihoods.

Author(s)

Dr. Lisa J Strug lisa.strug@utoronto.ca

References

Strug et. al, 2009. EJHG, 17:1171-1181;

See Also

evian_logit_plotsnp, li_plot, logit_snp, logit_plotsnp

Examples

## Not run: 
evian_logit(data, map, ycol, xcols, formula, model = "all", 
            rfdum = 'alphanum', m = 200, bse, lolim = log(0.025), 
            hilim = log(4), plot_li = TRUE, yaxislim = 10, 
            showmaxlr = 3, kcutoff = 32, grays = TRUE, 
            save_as_png = FALSE)
 
## End(Not run)

[Package evian version 1.1.0 Index]