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] 28.57616788 57.41942890 38.00290254 18.39142691 14.46354001 17.65108163
##   [7] 17.81531170 78.96611218 11.58485504 12.15293826  0.22816983 21.94881574
##  [13] 11.05548102 44.61881473 11.27345303 12.09194494  2.61502377 22.27845476
##  [19]  9.05528560 19.97219820  5.68075300  6.64137262 80.39038083 35.10028249
##  [25] 34.70130951 21.64564662 33.71500315 14.40791558  4.03888993 32.98578607
##  [31] 40.26523846 32.08567301 54.07060672 16.90044205 10.62234495 10.01908877
##  [37] 30.42880957 39.22976444 74.02452423 33.63132353 29.05642265  5.69187208
##  [43]  2.19404227  3.63666124 15.11180850 11.54513481  2.86048804 25.13240483
##  [49] 38.18396519 34.62223508 19.75676414 49.22415362 54.66114191  7.88081341
##  [55] 23.86850887 39.71984169 26.07085171  4.48332305 23.28883038 27.89417835
##  [61] 57.39949048 21.23334501  7.76252169 20.63809134  5.99230496 14.58871636
##  [67] 30.47650618  8.85502336  5.31155590 25.16787033 14.69519615  0.59886724
##  [73] 14.96898804 15.58692076 27.11110084 17.37507797  2.76784053 18.46605438
##  [79] 21.57233867 43.04767375 16.56548948 18.71833416  9.79845299 39.70599905
##  [85] 40.16866298 47.27440242  0.01559868 10.55773595  3.67692757 46.45958087
##  [91] 16.97636603  5.29303169 19.56677162 23.66629905  3.31333418 47.14609978
##  [97]  9.85571354 43.31992542 28.20963031 26.11180442