Compute Kernel Matrix

Description

Computes kernel matrix for a given genotype matrix

Usage

computeKernel(Z, kernel=c("linear.podkat", "quadratic.podkat",
              "localsim.podkat", "linear.SKAT", "quadratic.SKAT",
              "localsim.SKAT"), weights=NULL, pos=NULL, width=1000)

Arguments

Details

This function computes a kernel matrix for a given genotype matrix Z and a given kernel. It supposes that Z is a matrix-like object (a numeric matrix, a sparse matrix, or an object of class GenotypeMatrix) in which rows correspond to samples and columns correspond to variants. There are six different kernels available: “linear.podkat”, “quadratic.podkat”, “localsim.podkat”, “linear.SKAT”, “quadratic.SKAT”, and “localsim.SKAT”. All of these kernels can be used with or without weights. The weights can be specified with the weights argument which must be a numeric vector with as many elements as the matrix Z has columns. If no weighting should be used, weights must be set to NULL.

The position-dependent kernels “linear.podkat”, “quadratic.podkat”, and “localsim.podkat” require the positions of the variants in Z. So, if any of these three kernels is selected, the argument pos is mandatory and must be a numeric vector with as many elements as the matrix Z has columns.

If the pos argument is NULL and Z is a GenotypeMatrix object, the positions in variantInfo(Z) are taken. In this case, all variants need to reside on the same chromosome. If the variants in variantInfo(Z) are from multiple chromosomes, computeKernel quits with an error. As said, this only happens if pos is NULL, otherwise the pos argument has priority over the information stored in variantInfo(Z).

For details on how the kernels compute the pairwise similarities of genotypes, see Subsection 9.2 of the package vignette.

Value

a positive semi-definite kernel matrix with as many rows and columns as Z has rows

Author(s)

Ulrich Bodenhofer bodenhofer@bioinf.jku.at

References

http://www.bioinf.jku.at/software/podkat

Wu, M. C., Lee, S., Cai, T., Li, Y., Boehnke, M., and Lin, X. (2011) Rare-variant association testing for sequencing data with the sequence kernel association test. Am. J. Hum. Genet. 89, 82-93. DOI: 10.1016/j.ajhg.2011.05.029.

See Also

GenotypeMatrix

Examples

## create a toy example
A <- matrix(rbinom(50, 2, prob=0.2), 5, 10)
pos <- sort(sample(1:10000, ncol(A)))

## compute some unweighted kernels
computeKernel(A, kernel="linear.podkat", pos=pos, width=100)
computeKernel(A, kernel="localsim.podkat", pos=pos, width=100)
computeKernel(A, kernel="linear.SKAT")

## compute some weighted kernels
MAF <- colSums(A) / (2 * nrow(A))
weights <- betaWeights(MAF)
computeKernel(A, kernel="linear.podkat", pos=pos, weights=weights)
computeKernel(A, kernel="linear.SKAT", weights=weights)
computeKernel(A, kernel="localsim.SKAT", weights=weights)

Results

R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(podkat)
Loading required package: Rsamtools
Loading required package: GenomeInfoDb
Loading required package: stats4
Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:parallel':

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from 'package:stats':

    IQR, mad, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, append,
    as.data.frame, cbind, colnames, do.call, duplicated, eval, evalq,
    get, grep, grepl, intersect, is.unsorted, lapply, lengths, mapply,
    match, mget, order, paste, pmax, pmax.int, pmin, pmin.int, rank,
    rbind, rownames, sapply, setdiff, sort, table, tapply, union,
    unique, unsplit

Loading required package: S4Vectors

Attaching package: 'S4Vectors'

The following objects are masked from 'package:base':

    colMeans, colSums, expand.grid, rowMeans, rowSums

Loading required package: IRanges
Loading required package: GenomicRanges
Loading required package: Biostrings
Loading required package: XVector
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/podkat/computeKernel.Rd_%03d_medium.png", width=480, height=480)
> ### Name: computeKernel
> ### Title: Compute Kernel Matrix
> ### Aliases: computeKernel
> 
> ### ** Examples
> 
> ## create a toy example
> A <- matrix(rbinom(50, 2, prob=0.2), 5, 10)
> pos <- sort(sample(1:10000, ncol(A)))
> 
> ## compute some unweighted kernels
> computeKernel(A, kernel="linear.podkat", pos=pos, width=100)
     [,1] [,2] [,3] [,4] [,5]
[1,]    2 1.00 1.00    1    2
[2,]    1 9.36 2.36    2    3
[3,]    1 2.36 4.36    1    3
[4,]    1 2.00 1.00    2    3
[5,]    2 3.00 3.00    3    6
> computeKernel(A, kernel="localsim.podkat", pos=pos, width=100)
     [,1] [,2] [,3] [,4] [,5]
[1,] 1.00 0.62 0.77 0.90 0.80
[2,] 0.62 1.00 0.65 0.72 0.62
[3,] 0.77 0.65 1.00 0.77 0.77
[4,] 0.90 0.72 0.77 1.00 0.90
[5,] 0.80 0.62 0.77 0.90 1.00
> computeKernel(A, kernel="linear.SKAT")
     [,1] [,2] [,3] [,4] [,5]
[1,]    2    1    1    1    2
[2,]    1    9    2    2    3
[3,]    1    2    4    1    3
[4,]    1    2    1    2    3
[5,]    2    3    3    3    6
> 
> ## compute some weighted kernels
> MAF <- colSums(A) / (2 * nrow(A))
> weights <- betaWeights(MAF)
> computeKernel(A, kernel="linear.podkat", pos=pos, weights=weights)
             [,1]         [,2]         [,3]         [,4]         [,5]
[1,] 4.386439e+00 6.131825e-03 7.948981e-17 6.131804e-03 6.131804e-03
[2,] 6.131825e-03 1.115655e+01 1.145344e+00 2.532516e-05 1.761262e-01
[3,] 7.948981e-17 1.145344e+00 7.373779e+00 7.948981e-17 1.588755e-02
[4,] 6.131804e-03 2.532516e-05 7.948981e-17 2.530332e-05 2.530332e-05
[5,] 6.131804e-03 1.761262e-01 1.588755e-02 2.530332e-05 1.591286e-02
> computeKernel(A, kernel="linear.SKAT", weights=weights)
             [,1]         [,2]         [,3]         [,4]         [,5]
[1,] 3.976678e+00 4.951760e-17 4.951760e-17 4.951760e-17 9.903520e-17
[2,] 4.951760e-17 8.023070e+00 1.393797e-02 2.293961e-05 2.293961e-05
[3,] 4.951760e-17 1.393797e-02 4.004554e+00 4.951760e-17 1.393797e-02
[4,] 4.951760e-17 2.293961e-05 4.951760e-17 2.293961e-05 2.293961e-05
[5,] 9.903520e-17 2.293961e-05 1.393797e-02 2.293961e-05 1.396091e-02
> computeKernel(A, kernel="localsim.SKAT", weights=weights)
          [,1]      [,2]      [,3]      [,4]      [,5]
[1,] 1.0000000 0.9906608 0.9937739 0.9968978 0.9968869
[2,] 0.9906608 1.0000000 0.9906608 0.9937630 0.9937521
[3,] 0.9937739 0.9906608 1.0000000 0.9968761 0.9968869
[4,] 0.9968978 0.9937630 0.9968761 1.0000000 0.9999891
[5,] 0.9968869 0.9937521 0.9968869 0.9999891 1.0000000
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>