|
| 1 | +--- |
| 2 | +title: Using Eigen for eigenvalues |
| 3 | +author: Dirk Eddelbuettel |
| 4 | +license: GPL (>= 2) |
| 5 | +tags: eigen matrix |
| 6 | +summary: This example shows how to compute eigenvalues using Eigen |
| 7 | +layout: post |
| 8 | +src: 2013-01-11-eigen-eigenvalues.cpp |
| 9 | +--- |
| 10 | + |
| 11 | +[A previous post](../armadillo-eigenvalues) showed how to compute |
| 12 | +eigenvalues using the [Armadillo](http://arma.sf.net) library via RcppArmadillo. |
| 13 | + |
| 14 | +Here, we do the same using [Eigen](http://eigen.tuxfamily.org) and |
| 15 | +the RcppEigen package. |
| 16 | + |
| 17 | + |
| 18 | + |
| 19 | + |
| 20 | +{% highlight cpp %} |
| 21 | +#include <RcppEigen.h> |
| 22 | + |
| 23 | +// [[Rcpp::depends(RcppEigen)]] |
| 24 | + |
| 25 | +using Eigen::Map; // 'maps' rather than copies |
| 26 | +using Eigen::MatrixXd; // variable size matrix, double precision |
| 27 | +using Eigen::VectorXd; // variable size vector, double precision |
| 28 | +using Eigen::SelfAdjointEigenSolver; // one of the eigenvalue solvers |
| 29 | + |
| 30 | +// [[Rcpp::export]] |
| 31 | +VectorXd getEigenValues(Map<MatrixXd> M) { |
| 32 | + SelfAdjointEigenSolver<MatrixXd> es(M); |
| 33 | + return es.eigenvalues(); |
| 34 | +} |
| 35 | +{% endhighlight %} |
| 36 | + |
| 37 | + |
| 38 | +We can illustrate this easily via a random sample matrix. |
| 39 | + |
| 40 | +{% highlight r %} |
| 41 | +set.seed(42) |
| 42 | +X <- matrix(rnorm(4*4), 4, 4) |
| 43 | +Z <- X %*% t(X) |
| 44 | + |
| 45 | +getEigenValues(Z) |
| 46 | +{% endhighlight %} |
| 47 | + |
| 48 | + |
| 49 | + |
| 50 | +<pre class="output"> |
| 51 | +[1] 0.3319 1.6856 2.4099 14.2100 |
| 52 | +</pre> |
| 53 | + |
| 54 | + |
| 55 | +In comparison, R gets the same results (in reverse order) and also returns the eigenvectors. |
| 56 | + |
| 57 | +{% highlight r %} |
| 58 | +eigen(Z) |
| 59 | +{% endhighlight %} |
| 60 | + |
| 61 | + |
| 62 | + |
| 63 | +<pre class="output"> |
| 64 | +$values |
| 65 | +[1] 14.2100 2.4099 1.6856 0.3319 |
| 66 | + |
| 67 | +$vectors |
| 68 | + [,1] [,2] [,3] [,4] |
| 69 | +[1,] 0.69988 -0.55799 0.4458 -0.00627 |
| 70 | +[2,] -0.06833 -0.08433 0.0157 0.99397 |
| 71 | +[3,] 0.44100 -0.15334 -0.8838 0.03127 |
| 72 | +[4,] 0.55769 0.81118 0.1413 0.10493 |
| 73 | +</pre> |
| 74 | + |
| 75 | + |
| 76 | +Eigen has other _a lot_ of other decompositions, see [its documentation](http://eigen.tuxfamily.org/) |
| 77 | +for more details. |
0 commit comments