[Lme4-commits] r1477 - pkg/lme4Eigen/R

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Thu Dec 15 22:10:40 CET 2011


Author: dmbates
Date: 2011-12-15 22:10:39 +0100 (Thu, 15 Dec 2011)
New Revision: 1477

Added:
   pkg/lme4Eigen/R/GHrule.R
Log:
Added an R function, based on SparseGrid:::GQN to evaluate a matrix of nodes, weights and log-density for Gauss-Hermite quadrature against a Gaussian density.


Added: pkg/lme4Eigen/R/GHrule.R
===================================================================
--- pkg/lme4Eigen/R/GHrule.R	                        (rev 0)
+++ pkg/lme4Eigen/R/GHrule.R	2011-12-15 21:10:39 UTC (rev 1477)
@@ -0,0 +1,138 @@
+GQrule <- function (level) {
+    stopifnot(length(level) == 1,
+              (level <- as.integer(level)) > 0L,
+              level < 26L)
+    fr <- as.data.frame(switch(level,
+           list(n = 0, w = 1),
+           list(n = 1, w = 0.5),
+           list(n = c(0, 1.73205080756888), w = c(0.666666666666667, 0.166666666666667)),
+           list(n = c(0.741963784302726, 2.33441421833898),
+                w = c(0.454124145231931, 0.0458758547680685)),
+           list(n = c(0, 1.35562617997427, 2.85697001387281),
+                w = c(0.533333333333333, 0.222075922005613, 0.0112574113277207)),
+           list(n = c(0.616706590192594, 1.88917587775371, 3.32425743355212),
+                w = c(0.408828469556029, 0.0886157460419145, 0.00255578440205624)),
+           list(n = c(0, 1.15440539473997, 2.36675941073454, 3.75043971772574),
+                w = c(0.457142857142858, 0.240123178605013, 0.0307571239675865,
+                0.000548268855972219)),
+           list(n = c(0.539079811351375, 1.63651904243511, 2.80248586128754, 4.14454718612589),
+                w = c(0.373012257679077, 0.117239907661759, 0.00963522012078826,
+                0.000112614538375368)),
+           list(n = c(0, 1.02325566378913, 2.07684797867783, 3.20542900285647, 4.51274586339978),
+                w = c(0.406349206349207, 0.244097502894939, 0.049916406765218, 
+                0.00278914132123177, 2.23458440077466e-05)),
+           list(n = c(0.484935707515498, 1.46598909439116, 2.48432584163895, 
+                3.58182348355193, 4.85946282833231),
+                w = c(0.344642334932019, 0.135483702980267, 0.0191115805007703, 
+                0.00075807093431222, 4.31065263071831e-06)),
+           list(n = c(0, 0.928868997381064, 1.87603502015485, 2.86512316064364, 
+                3.93616660712998, 5.18800122437487),
+                w = c(0.36940836940837, 0.24224029987397, 0.0661387460710576, 
+                0.00672028523553727, 0.000195671930271223, 8.1218497902149e-07)),
+           list(n = c(0.444403001944139, 1.34037519715162, 2.2594644510008, 
+                3.2237098287701, 4.27182584793228, 5.50090170446775),
+                w = c(0.32166436151283, 0.14696704804533, 0.0291166879123641, 
+                0.00220338068753318, 4.83718492259061e-05, 1.49992716763716e-07)),
+           list(n = c(0, 0.85667949351945, 1.72541837958824, 2.62068997343221, 
+                3.56344438028163, 4.59139844893652, 5.8001672523865),
+                w = c(0.340992340992341, 0.237871522964136, 0.0791689558604501, 
+                0.0117705605059965, 0.000681236350442926, 1.15265965273339e-05,
+                2.7226276428059e-08)),
+           list(n = c(0.412590457954602, 1.24268895548546, 2.08834474570194, 
+                2.96303657983867, 3.88692457505977, 4.89693639734556, 
+                6.08740954690129),
+                w = c(0.302634626813019, 0.154083339842514, 0.0386501088242534, 
+                0.00442891910694741, 0.000200339553760744, 2.66099134406763e-06, 
+                4.86816125774839e-09)),
+           list(n = c(0, 0.799129068324548, 1.60671006902873, 2.43243682700976, 
+                3.28908242439877, 4.19620771126902, 5.19009359130478, 
+                6.36394788882984),
+                w = c(0.318259518259518, 0.232462293609732, 0.0894177953998444, 
+                0.0173657744921376, 0.00156735750354996, 5.64214640518902e-05, 
+                5.9754195979206e-07, 8.58964989963318e-10)),
+           list(n = c(0.386760604500557, 1.16382910055496, 1.95198034571633, 
+                2.7602450476307, 3.60087362417155, 4.49295530252001, 
+                5.47222570594934, 6.63087819839313),
+                w = c(0.286568521238012, 0.158338372750949, 0.0472847523540141, 
+                0.00726693760118474, 0.00052598492657391, 1.53000321624873e-05, 
+                1.30947321628682e-07, 1.49781472316183e-10)),
+           list(n = c(0, 0.751842600703896, 1.50988330779674, 2.28101944025299, 
+                3.07379717532819, 3.90006571719801, 4.77853158962998, 
+                5.74446007865941, 6.88912243989533),
+                w = c(0.299538370126608, 0.226706308468979, 0.0974063711627181, 
+                0.0230866570257112, 0.00285894606228465, 0.000168491431551339, 
+                4.01267944797987e-06, 2.80801611793058e-08, 2.58431491937492e-11)),
+           list(n = c(0.365245755507698, 1.0983955180915, 1.83977992150865, 
+                2.59583368891124, 3.37473653577809, 4.1880202316294, 
+                5.05407268544274, 6.0077459113596, 7.13946484914648),
+                w = c(0.272783234654288, 0.160685303893513, 0.0548966324802227, 
+                0.0105165177519414, 0.00106548479629165, 5.1798961441162e-05, 
+                1.02155239763698e-06, 5.90548847883655e-09, 4.41658876935871e-12)),
+           list(n = c(0, 0.71208504404238, 1.42887667607837, 2.15550276131694, 
+                2.89805127651575, 3.66441654745064, 4.46587262683103, 
+                5.32053637733604, 6.26289115651325, 7.38257902403043),
+                w = c(0.283773192751521, 0.220941712199144, 0.103603657276144, 
+                0.0286666910301185, 0.00450723542034204, 0.000378502109414268, 
+                1.53511459546667e-05, 2.53222003209287e-07, 1.22037084844748e-09, 
+                7.48283005405723e-13)),
+           list(n = c(0.346964157081356, 1.04294534880275, 1.74524732081413, 
+                2.45866361117237, 3.18901481655339, 3.94396735065732, 
+                4.73458133404606, 5.5787388058932, 6.51059015701366, 
+                7.61904854167976),
+                w = c(0.260793063449555, 0.161739333984, 0.061506372063976, 
+                0.013997837447101, 0.00183010313108049, 0.000128826279961929, 
+                4.40212109023086e-06, 6.12749025998296e-08, 2.48206236231518e-10, 
+                1.25780067243793e-13)),
+           list(n = c(0, 0.678045692440644, 1.35976582321123, 2.04910246825716, 
+                2.75059298105237, 3.46984669047538, 4.21434398168842, 
+                4.99496394478203, 5.82938200730447, 6.75144471871746, 
+                7.84938289511382),
+                w = c(0.270260183572877, 0.21533371569506, 0.108392285626419, 
+                0.0339527297865428, 0.00643969705140878, 0.000708047795481537, 
+                4.21923474255159e-05, 1.22535483614825e-06, 1.45066128449307e-08, 
+                4.97536860412175e-11, 2.09899121956567e-14)),
+           list(n = c(0.331179315715274, 0.995162422271216, 1.66412483911791, 
+                2.34175999628771, 3.03240422783168, 3.74149635026652, 
+                4.47636197731087, 5.24772443371443, 6.0730749511229, 
+                6.98598042401882, 8.07402998402171),
+                w = c(0.250243596586935, 0.161906293413675, 0.0671963114288899, 
+                0.0175690728808058, 0.00280876104757721, 0.000262283303255964, 
+                1.33459771268087e-05, 3.319853749814e-07, 3.36651415945821e-09, 
+                9.84137898234601e-12, 3.47946064787714e-15)),
+           list(n = c(0, 0.648471153534496, 1.29987646830398, 1.95732755293342, 
+                2.62432363405918, 3.30504002175297, 4.0047753217333, 
+                4.73072419745147, 5.49347398647179, 6.3103498544484, 
+                7.21465943505186, 8.29338602741735),
+                w = c(0.258509740808839, 0.209959669577543, 0.112073382602621, 
+                0.0388671837034809, 0.00857967839146566, 0.00116762863749786, 
+                9.3408186090313e-05, 4.08997724499215e-06, 8.77506248386172e-08, 
+                7.67088886239991e-10, 1.92293531156779e-12, 5.73238316780209e-16)),
+           list(n = c(0.317370096629452, 0.953421922932109, 1.59348042981642, 
+                2.24046785169175, 2.89772864322331, 3.56930676407356, 
+                4.26038360501991, 4.97804137463912, 5.7327471752512, 
+                6.54167500509863, 7.43789066602166, 8.50780351919526),
+                w = c(0.240870115546641, 0.161459512867, 0.0720693640171784, 
+                0.021126344408967, 0.00397660892918131, 0.000464718718779398, 
+                3.2095005652746e-05, 1.21765974544258e-06, 2.26746167348047e-08, 
+                1.71866492796487e-10, 3.71497415276242e-13, 9.39019368904192e-17)),
+           list(n = c(0, 0.622462279186076, 1.24731197561679, 1.87705836994784, 
+                2.51447330395221, 3.16277567938819, 3.82590056997249, 
+                4.50892992296729, 5.21884809364428, 5.9660146906067, 
+                6.76746496380972, 7.65603795539308, 8.71759767839959),
+                w = c(0.248169351176485, 0.20485102565034, 0.114880924303952, 
+                0.043379970167645, 0.0108567559914623, 0.0017578504052638, 
+                0.000177766906926527, 1.06721949052025e-05, 3.5301525602455e-07, 
+                5.73802386889938e-09, 3.79115000047719e-11, 7.10210303700393e-14, 
+                1.53003899799868e-17))
+           )
+    )
+    nr <- nrow(fr)
+    if (level %% 2L) {                 # level is odd
+        if (nr > 1L) {
+            fr <- rbind(fr[rev(2:nr), ], fr)
+        }
+    } else fr <- rbind(fr[rev(seq_len(nr)),], fr)
+    fr[seq_len(level %/% 2L), "n"] <- -fr[seq_len(level %/% 2L), "n"]
+    rownames(fr) <- NULL
+    as.matrix(within(fr, ldnorm <- dnorm(n, log=TRUE)))
+}



More information about the Lme4-commits mailing list