File: 6A-VectorMatrixAddon.Rd

package info (click to toggle)
fmultivar 240.10068-1
  • links: PTS
  • area: main
  • in suites: etch, etch-m68k
  • size: 1,492 kB
  • ctags: 272
  • sloc: fortran: 1,128; ansic: 764; sh: 22; makefile: 1
file content (603 lines) | stat: -rw-r--r-- 21,174 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
\name{VectorMatrixAddon}

\alias{VectorMatrixAddon}

\alias{triang}
\alias{Triang}
\alias{pascal}
\alias{colVec}
\alias{rowVec}
\alias{colIds}
\alias{rowIds}
\alias{colIds<-}
\alias{rowIds<-}

\alias{inv}
\alias{norm}
\alias{rk}
\alias{tr}
\alias{kron}
\alias{mexp}

\alias{vec}
\alias{vech}


\alias{tslag}
\alias{pdl}


\title{Vector/Matrix Arithmetics and Linear Algebra}


\description{
    
    A collection and description of functions 
    available for matrix arithmetics and linear 
    algebra. 
    \cr
    
    Only functions which have been added by Rmetrics 
    are documented. These functions are often useful 
    for the manipulation of multivariate time series. 
    \cr
    
    The origin of the functions is marked in the following way:
    
    R: part of R's base packages, \cr
    B: part of Rmetrics' fBasics package, \cr
    M: part of this Rmetrics package, fMultivar.
    \cr
    
    The functions are listed by topic.
    \cr
 
    General Matrix Functions:
    
    \tabular{ll}{  
    \code{matrix} \tab Creates a matrix from the given set of values, \cr
    \code{diag} \tab R: Creates a diagonal matrix or extracts diagonals, \cr
    \code{triang} \tab S: Extracs the lower tridiagonal part from a matrix, \cr
    \code{Triang} \tab S: Extracs the upper tridiagonal part from a matrix, \cr
    \code{pascal} \tab S: Creates a Pascal matrix, \cr
    \code{colVec} \tab S: Creates a column vector from a vector, \cr
    \code{rowVec} \tab S: Creates a row vector from a vector, \cr
    \code{as.matrix} \tab R: Attempts to turn its argument into a matrix, \cr
    \code{is.matrix} \tab R: Tests if its argument is a (strict) matrix, \cr
    \code{dimnames} \tab R: Retrieves or sets the dimnames of an object, \cr
    \code{colnames} \tab R: Retrieves or sets the column names, \cr 
    \code{colIds} \tab  R: ... use alternatively, \cr 
    \code{rownames} \tab  R: Retrieves or sets the row names, \cr   
    \code{rowIds} \tab  R: ... use alternatively. }
    
    Simple Matrix Operations:
            
    \tabular{ll}{ 
    \code{dim} \tab R: Returns the dimension of a matrix object, \cr
    \code{ncol} \tab R: Counts columns of a matrix object, \cr 
    \code{nrow} \tab R: Counts rows of a matrix object, \cr 
    \code{length} \tab R: Counts elements of a matrix object, \cr 
    \code{"["} \tab R: Subsets a matrix object, \cr
    \code{"[["} \tab R: Subsets a matrix object, \cr
    \code{cbind} \tab R: Augments a matrix object by columns, \cr 
    \code{rbind} \tab R: Augments a matrix object by rows. } 

    Basic Statistical Functions:
            
    \tabular{ll}{
    \code{var} \tab R: Returns the variance matrix, \cr
    \code{cov} \tab R: Returns the covariance matrix, \cr 
    \code{colStats} \tab B: Calculates column statistics, \cr
    \code{rowStats} \tab B: Calculates row statistics, \cr
    \code{colMeans} \tab B: Calculates column means, \cr 
    \code{rowMeans} \tab B: Calculates row means, \cr 
    \code{colAvgs} \tab B: Calculates column averages, \cr
    \code{rowAvgs} \tab B: Calculates row averages, \cr
    \code{colVars} \tab B: Calculates column variances, \cr 
    \code{rowVars} \tab B: Calculates row variances, \cr 
    \code{colStdevs} \tab B: Calculates column standard deviations, \cr
    \code{rowStdevs} \tab B: Calculates row standard deviations, \cr
    \code{colSkewness} \tab B: Calculates column skewness, \cr
    \code{rowSkewness} \tab B: Calculates row skewness, \cr
    \code{colKurtosis} \tab B: Calculates column kurtosis, \cr 
    \code{rowKurtosis} \tab B: Calculates row kurtosis, \cr 
    \code{colCumsums} \tab B: Calculates column cumulated sums, \cr  
    \code{rowCumsums} \tab B: Calculates row cumulated sums. }  
            
    Linear algebra:
    
    \tabular{ll}{
    \code{\%*\%} \tab R: Returns the product of two matrices, \cr
    \code{\%x\%}, \code{kron} \tab R: Returns the Kronecker product, \cr 
    \code{mexp} \tab M: Computes the exponential of a square matrix, \cr
    \code{det} \tab R: Returns the determinante of a matrix, \cr
    \code{inv} \tab S: Returns the inverse of a matrix, \cr
    \code{norm} \tab S: Returns the norm of a matrix, \cr 
    \code{rk} \tab S: Returns the rank of a matrix, \cr
    \code{tr} \tab S: Returns trace of a matrix, \cr 
    \code{t} \tab R: Returns the transposed matrix, \cr
    \code{vech} \tab M: Is the operator that stacks the lower triangle, \cr
    \code{vec} \tab M: Is the operator that stacks a matrix.}
    
    More linear algebra:
    
    \tabular{ll}{
    \code{chol} \tab R: Returns the Cholesky factor matrix, \cr
    \code{eigen} \tab R: Returns eigenvalues and eigenvectors, \cr 
    \code{svd} \tab R: Returns the singular value decomposition, \cr
    \code{kappa} \tab R: Returns the condition number of a matrix, \cr 
    \code{qr} \tab R: Returns the QR decomposition of a matrix, \cr
    \code{solve} \tab R: Solves a system of linear equations, \cr
    \code{backsolve} \tab R: ... use when the matrix is upper triangular, \cr
    \code{forwardsolve} \tab R: ... use when the matrix is lower triangular. } 
    
    Time Series Generation:
    
    \tabular{ll}{
    \code{tslag} \tab R: Lagged or leading vector/matrix of selected order(s), \cr
    \code{pdl} \tab R: Regressor matrix for polynomial distributed lags. } 
     

}


\usage{
triang(x)
Triang(x)
pascal(n)

colVec(x)
rowVec(x)
colIds(x, \dots)
rowIds(x, \dots)

inv(x)
norm(x, p = 2)
rk(x, method = c("qr", "chol"))
tr(x)
kron(x, y)

mexp(x, order = 8, method = c("pade", "taylor"))

vec(x)
vech(x)

pdl(x, d = 2, q = 3, trim = FALSE)
tslag(x, k = 1, trim = FALSE)
}


\arguments{

    \item{d}{
        [pdl] - \cr
        an integer specifying the order of the polynomial. 
        }
    \item{k}{
        [tslag] - \cr
        an integer value, the number of positions the new series is 
        to lag or to lead the input series. 
        }
    \item{method}{
        [mexp] - \cr
        Two methods are provided for the computation of the exponential
        of a matrix: Taylor series selected by \code{method="taylor"},
        and Pad\'e approximation with scaling and squaring to increase 
        precision selected by \code{method="pade"} which is the default
        method.
        \cr
        [rk] - \cr
        a character value, the dimension of the square matrix. One can
        choose from two methods: For \code{method = "qr"} the rank is
        computed as \code{qr(x)\$rank}, or alternatively for 
        \code{method="chol"} the rank is
        computed as \code{attr(chol(x, pivot=TRUE), "rank")}.
        }
    \item{n}{
        [pascal] - \cr
        an integer value, the dimension of the square matrix.
        }
    \item{order}{
        [mexp] - \cr
        The order of approximation to be used, an integer value.
        The best value for this depends on machine precision.
        }
    \item{p}{
        [norm] - \cr
        an integer value, \code{1}, \code{2} or \code{Inf}.
        \code{p=1} - The maximum absolute column sum norm which is defined 
        as the maximum of the sum of the absolute valued elements of columns 
        of the matrix.
        \code{p=2} - The spectral norm is "the norm" of a matrix \code{X}. 
        This value is computed as the square root of the maximum eigenvalue 
        of \code{CX} where \code{C} is the conjugate transpose.
        \code{p=Inf} - The maximum absolute row sum norm is defined 
        as the maximum of the sum of the absolute valued elements
        of rows of the matrix.
        }
    \item{q}{
        [pdl] - \cr
        an integer specifying the number of lags to use in 
        creating polynomial distributed lags. This must be 
        greater than d. 
        }
    \item{trim}{
        [pdl] - \cr
        a logical flag; if TRUE, the missing values at 
        the beginning of the returned matrix will be trimmed. 
        \cr
        [tslag] - \cr
         a logical flag, if TRUE, the missing values at the 
        beginning or end of the returned series will be trimmed. 
        The default value is FALSE. 
        }
    \item{x, y}{
        a numeric matrix.
        \cr
        [tslag] - \cr
        a numeric vector, missing values are allowed.
        \cr
        [pdl] - \cr
        a numeric vector.
        }
    \item{\dots}{
        arguments to be passed.
        }
        
}


\details{

    \bold{Function from R's Base Package:}
    \cr
    
    Most of the functions are described in their R help pages which we
    recommend to consult for further information. For the additiotnal 
    functions added by Rmetrics we give a brief introduction.
    \cr
    
    
    \bold{General Functions:}
    \cr
    
    Functions to generate matrices and related functions are described 
    in the help page \code{\link{matrix}}. To "decorate" these objects 
    several naming functions are available, a description can be found 
    on the help pages \code{\link{dimnames}} and \code{\link{rownames}}.
    \cr
    
    The function \code{pascal} generates a Pascal matrix of order \code{n} 
    which is a symmetric, positive, definite matrix with integer entries 
    made up from Pascal's triangle. The determinant of a Pascal matrix is 1. 
    The inverse of a Pascal matrix has integer entries. If \code{lambda} 
    is an eigenvalue of a Pascal matrix, then \code{1/lambda} is also an 
    eigenvalue of the matrix. Pascal matrices are ill-conditioned. 
    \cr
    
    The functions \code{triang} and \code{Triang} allow to transform a
    square matrix to a lower or upper triangular form. 
    A triangular matrix is either an upper triangular matrix or lower 
    triangular matrix. For the first case all matrix elements \code{a[i,j]}
    of matrix \code{A} are zero for \code{i>j}, whereas in the second case
    we have just the opposite situation. A lower triangular matrix is 
    sometimes also called left triangular. In fact, triangular matrices 
    are so useful that much computational linear algebra begins with 
    factoring or decomposing a general matrix or matrices into triangular 
    form. Some matrix factorization methods are the Cholesky factorization 
    and the LU-factorization. Even including the factorization step, 
    enough later operations are typically avoided to yield an overall 
    time savings. Triangular matrices have the following properties: the 
    inverse of a triangular matrix is a triangular matrix, the product of 
    two triangular matrices is a triangular matrix, the determinant of a 
    triangular matrix is the product of the diagonal elements, the 
    eigenvalues of a triangular matrix are the diagonal elements.
    \cr
    
    The functions \code{colVec} and \code{rowVec} transform a vector into 
    a column and row vector, respectively. A column vector is a matrix 
    object with one column, and a row vector is a matrix object with one row.
    \cr
    
    The functions \code{dimnames}, \code{colname}, and \code{rowname}
    can be used to retrieve and set the names of matrices. The functions 
    \code{rowIds}, \code{colIds}, are S-Plus like synonyms.
    \cr
    
    
    \bold{Simple Matrix Operations:}
    \cr
    
    The functions \code{\link{dim}}, \code{\link{nrow}} and 
    \code{\link{ncol}} are functions to extract the dimension and  
    the number of rows and columns of a matrix. 
    \cr
    
    The usual arithmetic operators, logical operators and mathematical 
    functions like \code{sqrt} or \code{exp} and \code{log}
    operate on matrices element by element. Note, that \code{"*"}
    is not matrix multiplication, instead we have to use \code{"\%*\%"}.
    \cr
    
    The methods \code{"["} and \code{"[["} are suited to extract
    subsets from a matrix, to delete rows and columns, or to permute
    rows and columns. 
    \cr
    
    
    \bold{Basic Statistical Functions:}
    \cr
    
    The functions \code{var} and \code{cov} compute the variance and
    covariance of a matrix.
    \cr
    
    Beside these functions \code{\link{colMeans}} and \code{\link{rowMeans}} 
    are R functions which compute the mean of columns and rows of a matrix. 
    Rmetrics has added further functions to compute column- or rowwise 
    variances, standard deviations, skewness, kurtosis and cumulated sums. 
    Two general  functions named \code{rowStats} and 
    \code{colStats} allow to apply through the argument list any 
    function to compute row and column statistics from matrices.
    \cr
    
    
    \code{Linear Algebra:}
    \cr
    
    Matrix multiplication is done using the operator \code{\%*\%}. 
    \cr
    
    The \emph{Kronecker product} can be computed using the operator 
    \code{\%x\%} or alternatively using the function \code{kron}.
    \cr
    
    The function \code{mexp} computes the exponential of a square 
    matrix \eqn{x}, defined as the sum from \eqn{r=0} to infinity of
    \eqn{x^r/r!}. Two methods are provided Taylor series and Pad\'e
    approximation with scaling and squaring to increase precision.
    Also reported, as the \code{"accuracy"} attribute of the result,
    is an upper bound for the L2 norm of the Cauchy error
    \code{mexp(a, order + 10, method) - mexp(a, order, method)}, and
    the used \code{method} and \code{order}. 

    The function \code{det} computes the determinant of a matrix.
    \cr
    
    The inverse of a square matrix \code{inv(X)} of dimension
    \code{n} is defined so that 
    \code{X \%*\% inv(X) = inv(X) \%*\% X = diag(n)} where
    the matrix \code{diag(n)} is the \code{n}-dimensional identity matrix.
    A precondition for the existence of the matrix inverse is that the
    determinant of the matrix \code{det(X)} is nonzero. 
    \cr

    The function \code{t} computes the transposed of a square matrix.
    \cr
    
    The function \code{vec} implements the operator that stacks a matrix
    as a column vector, to be more precise in a matrix with one column.
    vec(X) = (X11, X21, ..., XN1, X12, X22, ..., XNN).
    The function \code{vech} implements the operator that stacks the lower 
    triangle of a NxN matrix as an N(N+1)/2x1 vector:
    vech(X) =(X11, X21, X22, X31, ..., XNN), to be more precise in a 
    matrix with one row.
    \cr

    The function \code{tr} computes the trace of a square matrix which
    is the sum of the diagonal elements of the matrix under consideration.
    \cr
    
    The function \code{rk} computes the rank of a matrix which is 
    the dimension of the range of the matrix corresponding to the number 
    of linearly independent rows or columns of the matrix, or to the 
    number of nonzero singular values. The rank of a matrix is also
    named inear map. 
    \cr
    
    The function \code{norm} computes the norm of a matrix. Three choices 
    are possible: 
    \code{p=1} - The maximum absolute column sum norm which is defined 
    as the maximum of the sum of the absolute valued elements of columns 
    of the matrix.
    \code{p=2} - The spectral norm is "the norm" of a matrix \code{X}. 
    This value is computed as the square root of the maximum eigenvalue 
    of \code{CX} where \code{C} is the conjugate transpose.
    \code{p=Inf} - The maximum absolute row sum norm is defined 
    as the maximum of the sum of the absolute valued elements
    of rows of the matrix.
    \cr
    

    \code{More Linear Algebra:}
    \cr
    
    The function \code{chol} returns the Cholesky factor matrix,
    \code{eigen} returns eigenvalues and eigenvectors,
    \code{svd} returns the singular value decomposition,
    \code{kappa} estimate the condition number of a matrix,
    \code{qr} returns the QR decomposition of a matrix,
    \code{ginv} returns the Moore-Penrose generalized inverse,
    \code{solve} solves a system of linear equations, use
    \code{backsolve} when the matrix is upper triangular, and
    \code{forwardsolve} when the matrix is lower triangular.
    \cr
    
    
    \code{Time Series:}
    \cr
    The function \code{pdl} returns a regressor matrix suitable 
    for polynomial distributed lags. 
    \cr
    The function \code{tslag} returns a lagged/led vector or matrix 
    for given time series data.      
}


\references{

Higham N.J., (2002);
    \emph{Accuracy and Stability of Numerical Algorithms}, 
    2nd ed., SIAM.

Golub, van Loan, (1996);
    \emph{Matrix Computations}, 
    3rd edition. Johns Hopkins University Press. 
    
Ward, R.C., (1977);  
    \emph{Numerical computation of the matrix exponential with 
    accuracy estimate},
    SIAM J. Num. Anal. 14, 600--610.

Moler C., van Loan C., (2003);
    \emph{Nineteen dubious ways to compute the exponential of a matrix, 
    twenty-five years later},
    SIAM Review 45, 3--49.
}


\examples{
## SOURCE("fMultivar.6A-MatrixAddon")

## Create Pascal Matrix:
   P = pascal(3)
   P
   # Create lower triangle matrix
   L = triang(P)
   L
   # Extract diagonal part
   diag(P)
   
## Add/Subtract/Multiply/Divide:  
   X = P
   # Multiply matrix with a constant
   3 * X
   # Multiply two matrices elementwise
   X * P                     
   # Multiplies rows/columns of a matrix by a vector
   X \%*\% diag(P)            
   diag(P) \%*\% X           
       
## Operate on Subsets of a Matrix:
   n = 3; i = 2; j = 3
   D = diag(1:3)
   # Return the dimension of a matrix
   dim(P)                         
   # Get the last colum of a matrix
   P[, ncol(P)]                   
   # Delete a column of a matrix
   P[, -i]                      
   # Permute the columns of a matrix
   P[c(3, 1, 2), ]              
   # Augments matrix horizontally 
   cbind(P, D)                           
      
## Apply a function to all Elements of a Matrix: 
   # Return square root for each element
   sqrt(P)
   # Exponentiate the matrix elementwise
   exp(P)
   # Compute the median of each column
   apply(P, 2, "median") 
   # Test on all elements of a matrix       
   all( P > 2 )   
   # test on any element in a matrix                
   any( P > 2 )                  
     
## More Matrix Operations:
   # Return the product of two matrices
   P \%*\% D   
   # Return the Kronecker Product                     
   P \%x\% D                        
   # Return the transposed matrix
   t(P)                           
   # Return the inverse of a matrix
   inv(P)  
   # Return the norm of a matrix                      
   norm(P)                        
   # Return the determinante of a matrix
   det(P)                         
   # Return the rank of a matrix
   rk(P)                            
   # Return trace of a matrix
   tr(P)                          
   # Return the variance matrix
   var(P)     
   # Return the covariance matrix                   
   cov(P) 
   # Stack a matrix
   vec(P) 
   # Stack the lower triangle
   vech(P)
      
## Matrix Exponential:
   # Test case 1 from Ward (1977)
   test1 = t(matrix(c(
     4, 2, 0,
     1, 4, 1,
     1, 1, 4), 3, 3))
   mexp(test1)
   # Results on Power Mac G3 under Mac OS 10.2.8 
   #                    [,1]               [,2]               [,3]
   # [1,] 147.86662244637000 183.76513864636857  71.79703239999643
   # [2,] 127.78108552318250 183.76513864636877  91.88256932318409
   # [3,] 127.78108552318204 163.67960172318047 111.96810624637124
   # -- these agree with ward (1977, p608)
\dontrun{
   # A naive alternative to mexp, using spectral decomposition:
   mexp2 = function(matrix){
     z = eigen(matrix, sym = FALSE)
     Re(z$vectors \%*\% diag(exp(z$values)) \%*\% solve(z$vectors)) }
   mexp2(test1)
   # -- hopelessly inaccurate in all but the first column.
   # Test case 4 from Ward (1977)
   test4 <- structure(c(
     0, 0, 0, 0, 0, 0, 0, 0, 0, 1e-10,
     1, 0, 0, 0, 0, 0, 0, 0, 0, 0,
     0, 1, 0, 0, 0, 0, 0, 0, 0, 0,
     0, 0, 1, 0, 0, 0, 0, 0, 0, 0,
     0, 0, 0, 1, 0, 0, 0, 0, 0, 0,
     0, 0, 0, 0, 1, 0, 0, 0, 0, 0,
     0, 0, 0, 0, 0, 1, 0, 0, 0, 0,
     0, 0, 0, 0, 0, 0, 1, 0, 0, 0,
     0, 0, 0, 0, 0, 0, 0, 1, 0, 0,
     0, 0, 0, 0, 0, 0, 0, 0, 1, 0), .Dim = c(10, 10))
   attributes(mexp(test4))
   max(abs(mexp(test4) - mexp2(test4)))
   #[1] 8.746826694186494e-08
   # -- mexp2 is accurate to 7 d.p., whereas mexp to at least 14 d.p.              
}
   
## More Linear Algebra:
   X = P; b = c(1, 2, 3)
   # Return the Cholesky factor matrix
   chol(X)                        
   # Return eigenvalues and eigenvectors
   eigen(X)                       
   # Return the singular value decomposition
   svd(X)                         
   # Return the condition number of a matrix
   kappa(X)                       
   # Return the QR decomposition of a matrix
   qr(X)                          
   # Solve a system of linear equations
   # ... use backsolve when the matrix is upper triangular
   # ... use forwardsolve when the matrix is lower triangular
   solve(X, b)  
   backsolve(Triang(X), b)
   solve(Triang(X), b)                 
   forwardsolve(triang(X), b) 
   solve(triang(X), b)                                        
}


\author{

    Marina Shapira and David Firth for the \code{mexp} function, \cr
    Diethelm Wuertz for the Rmetrics \R-port.
    
}


\keyword{math}