linear_algebra_stats

Simple Linear Algebra Functions for Statistics

Description

‘broadcast’ provides some simple Linear Algebra Functions for Statistics:
cinv();
sd_lc().


Usage

cinv(x)

sd_lc(X, vc, bad_rp = NaN)

Arguments

x a real symmetric positive-definite square matrix.
X a numeric (or logical) matrix of multipliers/constants
vc the variance-covariance matrix for the (correlated) random variables.
bad_rp if vc is not a Positive (semi-) Definite matrix, give here the value to replace bad standard deviations with.

Details

cinv()
cinv() computes the Choleski inverse of a real symmetric positive-definite square matrix.

sd_lc()
Given the linear combination X %*% b, where:

  • X is a matrix of multipliers/constants;

  • b is a vector of (correlated) random variables;

  • vc is the symmetric variance-covariance matrix for b;

sd_lc(X, vc) computes the standard deviations for the linear combination X %*% b, without making needless copies.
sd_lc(X, vc) will use much less memory than a base ‘R’ approach.
sd_lc(X, vc) may possibly, but not necessarily, be faster than a base ‘R’ approach (depending on the Linear Algebra Library used for base ‘R’).


Value

For cinv():
A matrix.

For sd_lc():
A vector of standard deviations.

References

John A. Rice (2007), Mathematical Statistics and Data Analysis (6th Edition)

See Also

chol, chol2inv

Examples

library("broadcast")


vc <- datasets::ability.cov$cov
X <- matrix(rnorm(100), 100, ncol(vc))

solve(vc)
##              general       picture        blocks         maze       reading
## general  0.082259001 -0.0312406436 -7.750932e-03 -0.013309494 -2.061705e-02
## picture -0.031240644  0.2369906996 -2.484938e-02  0.017844845  8.603286e-04
## blocks  -0.007750932 -0.0248493822  1.344272e-02 -0.012544830 -3.802671e-05
## maze    -0.013309494  0.0178448450 -1.254483e-02  0.101625400  5.508423e-03
## reading -0.020617049  0.0008603286 -3.802671e-05  0.005508423  5.713620e-02
## vocab   -0.002420800  0.0019394999 -1.157864e-03 -0.002857265 -2.406969e-02
##                vocab
## general -0.002420800
## picture  0.001939500
## blocks  -0.001157864
## maze    -0.002857265
## reading -0.024069692
## vocab    0.020323179
cinv(vc) # faster than `solve()`, but only works on positive definite matrices
##              [,1]          [,2]          [,3]         [,4]          [,5]
## [1,]  0.082259001 -0.0312406436 -7.750932e-03 -0.013309494 -2.061705e-02
## [2,] -0.031240644  0.2369906996 -2.484938e-02  0.017844845  8.603286e-04
## [3,] -0.007750932 -0.0248493822  1.344272e-02 -0.012544830 -3.802671e-05
## [4,] -0.013309494  0.0178448450 -1.254483e-02  0.101625400  5.508423e-03
## [5,] -0.020617049  0.0008603286 -3.802671e-05  0.005508423  5.713620e-02
## [6,] -0.002420800  0.0019394999 -1.157864e-03 -0.002857265 -2.406969e-02
##              [,6]
## [1,] -0.002420800
## [2,]  0.001939500
## [3,] -0.001157864
## [4,] -0.002857265
## [5,] -0.024069692
## [6,]  0.020323179
all(round(solve(vc), 6) == round(cinv(vc), 6)) # they're the same
## [1] TRUE

sd_lc(X, vc)
##   [1] 43.4073171 41.2740401 53.3164093 14.0747582 35.2902031 88.2767914
##   [7] 24.8157933 12.0439739  7.1861024 28.6209981 40.2412223 41.0418633
##  [13] 38.6558326 12.6533232 82.3326336 57.2451051 26.3387797  1.7128853
##  [19]  2.3877058 41.7033251 18.5463942 41.1206122 78.3715210  6.5427717
##  [25] 17.5723607 27.9348138 27.1992269  5.3123103 37.3107653 15.7389135
##  [31]  6.5980134 57.3153458  0.7330632  4.7762495 62.9043613 20.7836610
##  [37] 45.5579373 20.9814090  3.5460098 45.2510950 28.2299643  7.5942212
##  [43] 13.3682944  1.1903421 13.8945040 14.3375901 19.0916791 30.6113497
##  [49] 48.2565760 15.6093831 12.5782772 15.9891274  6.9415761 38.5376621
##  [55]  8.9734555 16.5853854 74.5961937 21.8985844 41.4879770 50.5020260
##  [61] 66.1121175 15.0897261  6.1982754 50.5433333 25.5922920 12.3721057
##  [67] 32.5442976  6.5730409 16.2367350 27.3417614 30.8654273 20.8542828
##  [73] 10.0893138  2.8604279 25.8244753 30.8409941 50.0642030  3.8457165
##  [79] 33.8713509  4.7644235 23.4079907  7.6901473 17.2559830 21.0487438
##  [85] 30.8708882  4.9593600 28.6224839 22.6253557 27.7853181 13.7072525
##  [91]  0.3648497 16.6307393 10.5123965 15.5051164  0.1506426  3.8873052
##  [97] 54.0974694  5.5698790 10.3838945 37.9336070