expokit_dmexpv_Qmat: EXPOKIT dmexpv matrix exponentiation on Q matrix

Description Usage Arguments Details Value Author(s) See Also Examples

View source: R/rexpokit.R

Description

This function converts a matrix to COO format and exponentiates it via the EXPOKIT dmexpv function (designed for sparse matrices) and wrapper functions wrapalldmexpv_ around dmexpv.

Usage

1
2
3
4
  expokit_dmexpv_Qmat(Qmat = NULL, t = 2.1,
    inputprobs_for_fast = NULL, transpose_needed = TRUE,
    transform_to_coo_TF = TRUE, coo_n = NULL, anorm = NULL,
    check_for_0_rows = TRUE)

Arguments

Qmat

an input Q transition matrix

t

one or more time values to exponentiate by

inputprobs_for_fast

If NULL (default), the full probability matrix (Pmat) is returned. However, the full speed of EXPOKIT on sparse matrices will be exploited if inputprobs_for_fast=c(starting probabilities). In this case these starting probabilities are input to myDMEXPV directly, as v, and w, the output probabilities, are returned.

transpose_needed

If TRUE (default), matrix will be transposed (apparently EXPOKIT needs the input matrix to be transposed compared to normal)

transform_to_coo_TF

Should the matrix be tranposed to COO? COO format is required for EXPOKIT's sparse-matrix functions (like dmexpv and unlike the padm-related functions. Default TRUE; if FALSE, user must put a COO-formated matrix in Qmat. Supplying the coo matrix is probably faster for repeated calculations on large matrices.

coo_n

If a COO matrix is input, coo_n specified the order (# rows, equals # columns) of the matrix.

anorm

dmexpv requires an initial guess at the norm of the matrix. Using the

check_for_0_rows

If TRUE or a numeric value, the input Qmat is checked for all-zero rows, since these will crash the FORTRAN wrapalldmexpv function. A small nonzero value set to check_for_0_rows or the default (0.0000000000001) is input to off-diagonal cells in the row (and the diagonal value is normalized), which should fix the problem. R function norm might get slow with large matrices. If so, the user can input a guess manually (Lagrange seems to just use 1 or 0, if I recall correctly).

Details

From EXPOKIT:
* The method used is based on Krylov subspace projection
* techniques and the matrix under consideration interacts only
* via the external routine 'matvec' performing the matrix-vector
* product (matrix-free method).
*
* This is a customised version for Markov Chains. This means that a
* check is done within this code to ensure that the resulting vector
* w is a probability vector, i.e., w must have all its components
* in [0,1], with sum equal to 1. This check is done at some expense
* and the user may try DGEXPV which is cheaper since it ignores
* probability constraints.

COO (coordinated list) format is a compressed format that is
required for EXPOKIT's sparse-matrix functions (like dmexpv and
unlike EXPOKIT's padm-related functions.

COO (coordinated list) format is described here:

http://en.wikipedia.org/wiki/Sparse_matrix#Coordinate_list_.28COO.29

If Qmat is NULL (default), a default matrix is input.

Value

tmpoutmat the output matrix. wrapalldmexpv_ produces additional output relating to accuracy of the output matrix etc.; these can be by a direct call of dmexpv.

Author(s)

Nicholas J. Matzke [email protected] and Drew Schmidt [email protected]

See Also

mat2coo

expokit_wrapalldmexpv_tvals

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
# Example:
# Make a square instantaneous rate matrix (Q matrix)
# This matrix is taken from Peter Foster's (2001) "The Idiot's Guide
# to the Zen of Likelihood in a Nutshell in Seven Days for Dummies,
# Unleashed" at:
# \url{http://www.bioinf.org/molsys/data/idiots.pdf}
#
# The Q matrix includes the stationary base freqencies, which Pmat
# converges to as t becomes large.
Qmat = matrix(c(-1.218, 0.504, 0.336, 0.378, 0.126, -0.882, 0.252, 0.504, 0.168,
0.504, -1.05, 0.378, 0.126, 0.672, 0.252, -1.05), nrow=4, byrow=TRUE)

# Make a series of t values
tvals = c(0.001, 0.005, 0.01, 0.05, 0.1, 0.5, 1, 2, 5, 14)

# Exponentiate each with EXPOKIT's dmexpv (should be fast for large sparse matrices)
for (t in tvals)
	{
	Pmat = expokit_dmexpv_Qmat(Qmat=Qmat, t=t, transpose_needed=TRUE)
	cat("\n\nTime=", t, "\n", sep="")
	print(Pmat)
	}

rexpokit documentation built on Aug. 16, 2017, 1:04 a.m.