| Title: | Variance Heterogeneity Genome-wide Association Study - Reimplementation |
|---|---|
| Description: | The package reimplementation provides models and tests for variance heterogeneity genome-wide association study (vGWAS). |
| Authors: | Xia Shen |
| Maintainer: | Kristian Ullrich <[email protected]> |
| License: | GPL |
| Version: | 2025.05.18 |
| Built: | 2026-05-13 06:41:11 UTC |
| Source: | https://github.com/kullrich/vgwas |
The function performs the robust Brown-Forsythe test using the group medians.
bfmedian.test(formula, data, alpha = 0.05, na.rm = TRUE, verbose = TRUE)bfmedian.test(formula, data, alpha = 0.05, na.rm = TRUE, verbose = TRUE)
formula |
a formula of the form lhs ~ rhs where lhs gives the sample values and rhs the corresponding groups. |
data |
a tibble or data frame containing the variables in formula. |
alpha |
the level of significance to assess the statistical difference. Default is set to alpha = 0.05. |
na.rm |
a logical value indicating whether NA values should be stripped before the computation proceeds. Default us set to TRUE. |
verbose |
a logical for printing output to R console. |
Levene (1960) proposed a test for homogeneity of variances in k groups which is based on the ANOVA statistic applied to absolute deviations of observations from the corresponding group mean. The robust Brown-Forsythe version of the Levene-type test substitutes the group mean by the group median in the classical Levene statistic.
A list with class "owt" containing the following components:
statistic |
the Brown-Forsythe test statistic. |
parameter |
the parameter(s) of the approximate F distribution of the test statistic. |
p.value |
the p-value of the test. |
alpha |
the level of significance to assess the statistical difference. |
method |
the character string "Brown-Forsythe-Median Test". |
data |
a data frame containing the variables in which NA values (if exist) are removed. |
formula |
a formula of the form lhs ~ rhs where lhs gives the sample values and rhs the corresponding groups. |
Modified from the onewaytests package and vGWAS.
Kristian Ullrich
Brown, M. B. and Forsythe, A.B. (1974). Robust tests for
equality of variances.
Journal of the American Statistical Association, 69,
364-367.
Levene, H. (1960).
Robust Tests for Equality of Variances, in Contributions
to Probability and Statistics, ed. I. Olkin, Palo Alto, CA: Stanford Univ.
Press.
data(pheno) data(geno) df <- data.frame(phenotype = pheno, genotype = as.factor(geno[, 911])) bfmedian.test(phenotype ~ genotype, data = df)data(pheno) data(geno) df <- data.frame(phenotype = pheno, genotype = as.factor(geno[, 911])) bfmedian.test(phenotype ~ genotype, data = df)
The function performs the robust Brown-Forsythe test using the group medians. Instead of the ANOVA statistic, the Kruskal-Wallis ANOVA may also be applied using this function.
brown.forsythe.test(y, group, kruskal.test=FALSE)brown.forsythe.test(y, group, kruskal.test=FALSE)
y |
a numeric vector of data values. |
group |
factor of the data. |
kruskal.test |
a |
Levene (1960) proposed a test for homogeneity of variances in k groups which is based on the ANOVA statistic applied to absolute deviations of observations from the corresponding group mean. The robust Brown-Forsythe version of the Levene-type test substites the group mean by the group median in the classical Levene statistic.
A list with the following numeric components.
statistic |
the value of the test statistic. |
p.value |
the p-value of the test. |
method |
type of test performed. |
data.name |
a character string giving the name of the data. |
Modified from the lawstat package.
Xia Shen
Brown, M. B. and Forsythe, A.B. (1974). Robust tests for
equality of variances.
Journal of the American Statistical Association,
69, 364-367.
Levene, H. (1960).
Robust Tests for Equality of Variances, in Contributions
to Probability and Statistics, ed. I. Olkin, Palo Alto, CA: Stanford Univ.
Press.
data(pheno) data(geno) brown.forsythe.test(pheno, geno[,911])data(pheno) data(geno) brown.forsythe.test(pheno, geno[,911])
Chromosome indices for the markers of the simulated data
data(chr)data(chr)
A numeric vector of chromosome indices for the 20K simulated markers.
data("chr", package="vGWAS") table(chr)data("chr", package="vGWAS") table(chr)
The marker genotypes of the simulated data
data(geno)data(geno)
A character matrix of size (number of individuals) times (number of markers in the genome).
Note that there is only one column for each marker.
data("geno", package="vGWAS")data("geno", package="vGWAS")
The marker genotypes of the simulated data
data(geno.df)data(geno.df)
A data.frame of size (number of markers in the genome) times (number of individuals).
Note that there is only one column for each marker.
data("geno.df", package="vGWAS")data("geno.df", package="vGWAS")
The marker genotypes of the simulated data
data(geno.num)data(geno.num)
A matrix of size (number of individuals) times (number of markers in the genome).
Note that there is only one column for each marker.
data("geno.num", package="vGWAS")data("geno.num", package="vGWAS")
The marker genotypes of the simulated data
data(geno.sparse)data(geno.sparse)
A sparse matrix of size (number of markers in the genome) times (number of individuals).
Note that there is only one column for each marker.
data("geno.sparse", package="vGWAS")data("geno.sparse", package="vGWAS")
Calculates minor-allele-frequency (MAF)
getMAF(geno.matrix, geno.snp = "row", include.het = FALSE)getMAF(geno.matrix, geno.snp = "row", include.het = FALSE)
geno.matrix |
a |
geno.snp |
if individuals at columns and markers at rows use "row" else if individuals at rows and markers at columns use "col" |
include.het |
specify if heterozygous calls should be split and added equally to homozygous ref and alt counts (default = FALSE) |
a vector containing minor-allele-frequency values
Kristian Ullrich
data(geno) maf <- getMAF(geno.matrix = geno, geno.snp = "col") data(geno.num) maf.num <- getMAF(geno.matrix = geno.num, geno.snp = "col") data(geno.df) maf.df <- getMAF(geno.matrix = geno.df, geno.snp = "row") data(geno.sparse) maf.sparse <- getMAF(geno.matrix = geno.sparse, geno.snp = "row")data(geno) maf <- getMAF(geno.matrix = geno, geno.snp = "col") data(geno.num) maf.num <- getMAF(geno.matrix = geno.num, geno.snp = "col") data(geno.df) maf.df <- getMAF(geno.matrix = geno.df, geno.snp = "row") data(geno.sparse) maf.sparse <- getMAF(geno.matrix = geno.sparse, geno.snp = "row")
Map positions for the markers of the simulated data
data(map)data(map)
A numeric vector of chromosomal map positions of the 20K simulated markers.
data("map", package="vGWAS")data("map", package="vGWAS")
The package provides models and tests for variance heterogeneity genome-wide association study (vGWAS). "_PACKAGE"
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
Ronnegard, L., Shen, X. and Alam, M. (2011):
hglm: A Package for Fitting Hierarchical Generalized
Linear Models. The R Journal, 2(2), 20-28.
Brown, M. B. and Forsythe, A.B. (1974).
Robust tests for equality of variances. Journal
of the American Statistical Association, 69, 364-367.
Levene, H. (1960). Robust Tests for Equality
of Variances, in Contributions to Probability and Statistics,
ed. I. Olkin, Palo Alto, CA: Stanford Univ. Press.
R package lawstat for other types of nonparametric
variance tests and onewaystats.
Phenotypic values for the markers of the simulated data
data(pheno)data(pheno)
A numeric vector of the phenotypic values of 93 simulated individuals.
Note that there is only one column for each marker.
data("pheno", package="vGWAS") hist(pheno, breaks = 30)data("pheno", package="vGWAS") hist(pheno, breaks = 30)
The function plots the variance GWA result for the given scan object.
## S3 method for class 'vGWAS' plot(x, sig.threshold = NULL, low.log.p = 0, pch = 16, cex = 0.6, col.manhattan = c("slateblue4", "olivedrab"), col.sig.threshold = "darkgoldenrod", ...)## S3 method for class 'vGWAS' plot(x, sig.threshold = NULL, low.log.p = 0, pch = 16, cex = 0.6, col.manhattan = c("slateblue4", "olivedrab"), col.sig.threshold = "darkgoldenrod", ...)
x |
a result object from |
sig.threshold |
a numeric value giving the significance
threshold for |
low.log.p |
a numeric value giving the lower limit of the
|
pch |
point character. See |
cex |
size of points. See |
col.manhattan |
two colors as a vector for the Manhattan plot. |
col.sig.threshold |
one color for the significance threshold |
... |
not in use |
a plot for viewing vGWAS result.
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)
The function summarized the variance GWA result for the given scan object.
## S3 method for class 'vGWAS' summary(object, nrMarkers = 10, ...)## S3 method for class 'vGWAS' summary(object, nrMarkers = 10, ...)
object |
a result object from |
nrMarkers |
a numeric value giving the number of top markers to be summarized. |
... |
not in use |
a summary for viewing vGWAS result.
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
Ronnegard, L., Shen, X. and Alam, M. (2010):
hglm: A Package for Fitting Hierarchical Generalized
Linear Models. The R Journal, 2(2), 20-28.
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)
Variance Genome-wide association for using nonparametric variance test
vGWAS(phenotype, geno.matrix, kruskal.test = FALSE, marker.map = NULL, chr.index = NULL, pB = TRUE)vGWAS(phenotype, geno.matrix, kruskal.test = FALSE, marker.map = NULL, chr.index = NULL, pB = TRUE)
phenotype |
a |
geno.matrix |
a |
kruskal.test |
a |
marker.map |
a |
chr.index |
a |
pB |
show progress bar |
a data.frame containing columns of marker names,
chromosome indices, marker.map positions,
test statistic values, and p.value for each position.
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
Ronnegard, L., Shen, X. and Alam, M. (2010):
hglm: A Package for Fitting Hierarchical Generalized
Linear Models. The R Journal, 2(2), 20-28.
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)
The function does genomic control for the variance GWA result object.
vGWAS.gc(object, plot = TRUE, proportion = 1, ...)vGWAS.gc(object, plot = TRUE, proportion = 1, ...)
object |
a result object from |
plot |
a logical value turning on/off the QQ plot for genomic control. |
proportion |
a numeric value between 0 and 1 giving the proportion of obtained p-values to be used for genomic control. |
... |
not in use |
lambda |
estimated inflation ratio. |
lambda.se |
standard error of the estimated inflation ratio. |
gc.p.value |
p-values after genomic control. |
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])]) # ----- genomic control ----- # vgwa2 <- vGWAS.gc(vgwa) plot(vgwa2) summary(vgwa2)
The function calculates and reports the variance explained for a single marker by fitting a double generalized linear model. It gives both the variance explained by the mean and variance parts of model.
vGWAS.variance(phenotype, marker.genotype, print = TRUE)vGWAS.variance(phenotype, marker.genotype, print = TRUE)
phenotype |
a |
marker.genotype |
a |
print |
a |
The Value will only be available if
only.print = FALSE.
variance.mean |
the variance explained by the mean part of model. |
variance.disp |
the variance explained by the variance part of model. |
Xia Shen
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
package-vGWAS,
vGWAS, plot.vGWAS
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])])# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWAS(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, pB = FALSE) # ----- visualize the scan ----- # plot(vgwa) summary(vgwa) # ----- calculate the variance explained by the strongest marker ----- # vGWAS.variance(phenotype = pheno, marker.genotype = geno[, vgwa[["p.value"]] == min(vgwa[["p.value"]])])
Variance Genome-wide association for using nonparametric variance test and other
vGWASparallel(phenotype, geno.matrix, marker.map = NULL, chr.index = NULL, geno.snp = "row", method = "bfmedian", test.alpha = 0.05, test.na.rm = TRUE, p.adjust.method = "none", include.het = FALSE, pB = TRUE, ncores = 1)vGWASparallel(phenotype, geno.matrix, marker.map = NULL, chr.index = NULL, geno.snp = "row", method = "bfmedian", test.alpha = 0.05, test.na.rm = TRUE, p.adjust.method = "none", include.het = FALSE, pB = TRUE, ncores = 1)
phenotype |
a |
geno.matrix |
a |
marker.map |
a |
chr.index |
a |
geno.snp |
if individuals at columns and markers at rows use "row" else if individuals at rows and markers at columns use "col" |
method |
the test method to use (default = bfmedian). Default is set to the Brown-Forsythe's Test of Equality of Variances using group medians. There are 31 other tests available via the onewaytests package: Alvandi's F test ("af"), Alexander-Govern test ("ag"), Alvandi's generalized p-value ("agp"), One-way analysis of variance ("aov"), Approximate F test ("ap"), Adjusted Welch's heteroscedastic F test ("aw"), B square test ("b2"), Brown-Forsythe test ("bf"), Box F test ("box"), Cochran test ("cochran"), Generalized tests equivalent to Parametric Bootstrap ("gtb"), Generalized tests equivalent to Fiducial tests ("gtf"), Variance homogeneity tests ("homog"), James second order test ("james"), Johansen F test ("johansen"), Kruskal-Wallis test ("kw"), Modified Brown-Forsythe test ("mbf"), Mann-Whitney U test ("mw"), Anderson-Darling normility test ("nor_ad"), Cramer-vin Mises normility test ("nor_cvm"), Kolmogorov-Smirnov normility test ("nor_ks"), Pearson Chi-square normility test ("nor_pct"), Shapiro-Wilk normility test ("nor_sw"), Shapiro-Francia normility test ("nor_sf"), Permutation F test ("pf"), Scott-Smith test ("ss"), Student's t-test ("st"), Welch-Aspin test ("wa"), Welch's heteroscedastic F test with trimmed means and Winsorized variances ("welch"), Weerahandi's generalized F test ("wgf"), Welch's t-test ("wt"). |
test.alpha |
the level of significance to assess the statistical difference. Default is set to alpha = 0.05. |
test.na.rm |
a logical value indicating whether NA values should be stripped before the computation proceeds. Default us set to TRUE. |
p.adjust.method |
correction method (default = "none"). There are 8 p-value correction methods available via the p.adjust function: "holm", "hochberg", "hommel", "bonferroni", "BH", "BY", "fdr", "none" |
include.het |
specify if heterozygous calls should be split and added equally to homozygous ref and alt counts (default = FALSE) |
pB |
show progress bar |
ncores |
number of cores to parallelize (default = 1) |
a data.frame containing columns of marker names,
chromosome indices, marker.map positions,
test statistic values, and p.value for each position.
Xia Shen
Kristian Ullrich
Shen, X., Pettersson, M., Ronnegard, L. and Carlborg, O.
(2011): Inheritance beyond plain heritability:
variance-controlling genes in Arabidopsis thaliana.
PLoS Genetics, 8, e1002839.
Ronnegard, L., Shen, X. and Alam, M. (2010):
hglm: A Package for Fitting Hierarchical Generalized
Linear Models. The R Journal, 2(2), 20-28.
package-vGWAS onewaytests
# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", pB = FALSE) # ----- other test GWA scan ----- # vgwa.mw <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", method = "mw", pB = FALSE) # ----- multiple cores ----- # vgwa.st <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", method = "st", ncores = 2, pB = FALSE)# ----- load data ----- # data(pheno) data(geno) data(chr) data(map) # ----- variance GWA scan ----- # vgwa <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", pB = FALSE) # ----- other test GWA scan ----- # vgwa.mw <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", method = "mw", pB = FALSE) # ----- multiple cores ----- # vgwa.st <- vGWASparallel(phenotype = pheno, geno.matrix = geno, marker.map = map, chr.index = chr, geno.snp = "col", method = "st", ncores = 2, pB = FALSE)