netsnp {netgwas} | R Documentation |
Reconstructs intra- and inter- chromosomal conditional interactions among genetic loci
Description
This is one of the main functions of the netgwas package. This function can be used to reconstruct the intra- and inter-chromosomal interactions among genetic loci in diploids and polyploids. The input data can be belong to any biparental genotype data which contains at least two genotype states. Two methods are available to reconstruct the network, namely (1) approximation method, and (2) gibbs sampling within the Gaussian copula graphical model. Both methods are able to deal with missing genotypes.
Usage
netsnp(data, method = "gibbs", rho = NULL, n.rho = NULL, rho.ratio = NULL,
ncores = 1, em.iter = 5, em.tol = .001, verbose = TRUE)
Arguments
data |
An ( |
method |
Reconstructs intra- and inter- chromosomal conditional interactions (epistatic selection) network with three methods: "gibbs", "approx", and "npn". For a medium (~500) and a large number of variables we would recommend to choose "gibbs" and "approx", respectively. For a very large number of variables (> 2000) choose "npn". The default method is "gibbs". |
rho |
A decreasing sequence of non-negative numbers that control the sparsity level. Leaving the input as |
n.rho |
The number of regularization parameters. The default value is |
rho.ratio |
Determines the distance between the elements of |
ncores |
The number of cores to use for the calculations. Using |
em.iter |
The number of EM iterations. The default value is 5. |
em.tol |
A criteria to stop the EM iterations. The default value is .001. |
verbose |
Providing a detail message for tracing output. The default value is |
Details
Viability is a phenotype that can be considered. This function detects the conditional dependent short- and long-range linkage disequilibrium structure of genomes and thus reveals aberrant marker-marker associations that are due to epistatic selection. This function can be used to estimate conditional independence relationships between partially observed data that not follow Gaussianity assumption (e.g. continuous non-Gaussian, discrete, or mixed dataset).
Value
An object with S3 class "netgwas"
is returned:
Theta |
A list of estimated p by p precision matrices that show the conditional independence relationships patterns among genetic loci. |
path |
A list of estimated p by p adjacency matrices. This is the graph path corresponding to |
ES |
A list of estimated p by p conditional expectation corresponding to |
Z |
A list of n by p transformed data based on Gaussian copula. |
rho |
A |
loglik |
A |
data |
The |
Note
This function estimates a graph path . To select an optimal graph please refer to selectnet
.
Author(s)
Pariya Behrouzi and Ernst C. Wit
Maintainers: Pariya Behrouzi pariya.behrouzi@gmail.com
References
1. Behrouzi, P., and Wit, E. C. (2019). Detecting epistatic selection with partially observed genotype data by using copula graphical models. Journal of the Royal Statistical Society: Series C (Applied Statistics), 68(1), 141-160.
2. Behrouzi, P., and Wit, E. C. (2017c). netgwas: An R Package for Network-Based Genome-Wide Association Studies. arXiv preprint, arXiv:1710.01236.
3. D. Witten and J. Friedman. New insights and faster computations for the graphical lasso. Journal of Computational and Graphical Statistics, to appear, 2011.
4. Guo, Jian, et al. "Graphical models for ordinal data." Journal of Computational and Graphical Statistics 24.1 (2015): 183-204.
See Also
Examples
data(CviCol)
out <- netsnp(CviCol); out
plot(out)
#select optimal graph
epi <- selectnet(out)
plot(epi, vis="CI", xlab="markers", ylab="markers",
n.mem = c(24,14,17,16,19), vertex.size=4)
#Visualize interactive plot of the selected network
#Different colors for each chromosome
cl <- c(rep("red", 24), rep("white",14), rep("tan1",17),
rep("gray",16), rep("lightblue2",19))
plot(epi, vis="interactive", vertex.color= cl)
#Partial correlations between markers on genome
image(as.matrix(epi$par.cor), xlab="markers", ylab="markers", sub="")