[Rcpp-commits] r3686 - pkg/RcppEigen/inst/doc

[noreply at r-forge.r-project.org]

Author: dmbates
Date: 2012-07-10 20:27:15 +0200 (Tue, 10 Jul 2012)
New Revision: 3686

Modified:
   pkg/RcppEigen/inst/doc/code.R
Log:
Slight update in code (not actually included in vignette).


Modified: pkg/RcppEigen/inst/doc/code.R
===================================================================
--- pkg/RcppEigen/inst/doc/code.R	2012-07-10 18:26:14 UTC (rev 3685)
+++ pkg/RcppEigen/inst/doc/code.R	2012-07-10 18:27:15 UTC (rev 3686)
@@ -14,11 +14,16 @@
 typedef Map<MatrixXd>  MapMatd;
 typedef Map<MatrixXi>  MapMati;
 typedef Map<VectorXd>  MapVecd;
-inline MatrixXd AtA(const MapMatd& A) {
+inline  MatrixXd AtA(const MatrixXd& A) {
     int    n(A.cols());
     return   MatrixXd(n,n).setZero().selfadjointView<Lower>()
              .rankUpdate(A.adjoint());
 }
+inline  MatrixXd AAt(const MatrixXd& A) {
+    int    n(A.cols());
+    return   MatrixXd(n,n).setZero().selfadjointView<Lower>()
+             .rankUpdate(A);
+}
 '
 
 
@@ -123,16 +128,15 @@
 const VectorXd   fitted(X * betahat);
 const VectorXd    resid(y - fitted);
 const int            df(n - p);
-const double          s(resid.norm() / std::sqrt(double(df)));
-const VectorXd       se(s * llt.matrixL().solve(MatrixXd::Identity(p, p))
-                        .colwise().norm());
+const double        ssq(resid.squaredNorm() / double(df));
+const MatrixXd     vcov(ssq * llt.solve(MatrixXd::Identity(p, p)));
 return     List::create(Named("coefficients")   = betahat,
                         Named("fitted.values")  = fitted,
                         Named("residuals")      = resid,
-                        Named("s")              = s,
+                        Named("s")              = sqrt(ssq),
                         Named("df.residual")    = df,
                         Named("rank")           = p,
-                        Named("Std. Error")     = se);
+                        Named("vcov")           = vcov);
 '
 
 lltLS <- cxxfunction(signature(XX = "matrix", yy = "numeric"),
@@ -142,10 +146,9 @@
 str(lmFit <- with(trees, lm.fit(cbind(1, log(Girth)), log(Volume))))
 for (nm in c("coefficients", "residuals", "fitted.values", "rank"))
     stopifnot(all.equal(lltFit[[nm]], unname(lmFit[[nm]])))
-stopifnot(all.equal(lltFit[["Std. Error"]],
-                    unname(coef(summary(lm(log(Volume) ~ log(Girth), trees)))[,2])))
+stopifnot(all.equal(unname(vcov(lm(log(Volume) ~ log(Girth), trees))),
+                    lltVCFit$vcov))
 
-
 ## section 4.3
 
 dd <- data.frame(f1 = gl(4, 6, labels = LETTERS[1:4]),