[Splm-commits] r202 - pkg/R

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

Author: the_sculler
Date: 2015-11-29 23:04:31 +0100 (Sun, 29 Nov 2015)
New Revision: 202

Modified:
   pkg/R/sphtest.R
Log:
Some fixes to sphtest, still suboptimal.



Modified: pkg/R/sphtest.R
===================================================================
--- pkg/R/sphtest.R	2015-11-29 12:22:34 UTC (rev 201)
+++ pkg/R/sphtest.R	2015-11-29 22:04:31 UTC (rev 202)
@@ -6,8 +6,8 @@
 sphtest.formula <- function (x, data, index = NULL, listw,
                              spatial.model = c("lag", "error", "sarar"),
                              method = c("ML", "GM"), errors = c("KKP", "BSK"),...) {
-    ## performs a Hausman test of a FE model with spatial lag or error
-    ## against "alternative" with same spatial specification
+    ## performs a Hausman test of a RE model with spatial lag or error
+    ## against FE "alternative" with same spatial specification
 
     switch(match.arg(spatial.model),
     lag = {
@@ -29,10 +29,16 @@
 
     method <- switch(match.arg(method),
                      ML = {
+                         ## adapt argument
+                         spatial.error <- if(spatial.error) {
+                             spatial.error <- if(errors=="BSK") "b" else "kkp"
+                         } else {
+                             spatial.error <- "none"
+                         }
                          femod <- spml(x, data = data, index = index, listw = listw, lag = lag,
-                                       spatial.error = spatial.error, model = "within", errors = errors)
+                                       spatial.error = spatial.error, model = "within")
                          remod <- spml(x, data = data, index = index, listw = listw, lag = lag,
-                                       spatial.error = spatial.error, model = "random", errors = errors)
+                                       spatial.error = spatial.error, model = "random")
                          },
                      GM = {
                          femod <- spgm(x, data = data, index = index, listw = listw, lag = lag,
@@ -46,29 +52,60 @@
     }
 
 sphtest.splm <- function (x, x2, ...){
-    ## check that the models have the same specification but different effects
-    if (!all.equal(x$legacy, x2$legacy)) stop("The models are different")
-    if(x$ef.sph == x2$ef.sph) stop("Effects should be different")
+    ## check whether the models have been estimated by GM (different slots...)
+    is.gm <- !is.null(x$ef.sph)
 
-    ran <- match("random", c(x$ef.sph, x2$ef.sph))
-    if(ran == 1){
-	xwith <- x2
-	xbetw <- x
-	}
-    if(ran == 2){
-	xwith <- x
-	xbetw <- x2
-	}
+    if(is.gm) {
+        ## check that the models have the same specification but different effects
+        if (!all.equal(x$legacy, x2$legacy)) stop("The models are different")
+        if(x$ef.sph == x2$ef.sph) stop("Effects should be different")
 
-    ## test on coefficients (excluding SAR)
-    ## model order is irrelevant
+        ran <- match("random", c(x$ef.sph, x2$ef.sph))
+        if(ran == 1){
+            xwith <- x2
+            xbetw <- x
+        }
+        if(ran == 2){
+	    xwith <- x
+	    xbetw <- x2
+        }
+    
+      ## test on coefficients (excluding SAR)
+      ## model order is irrelevant
   
-    tc <- match(names(coef(xwith)), names(coef(xbetw)) )
+      tc <- match(names(coef(xwith)), names(coef(xbetw)) )
 
-    coef.wi <- coef(xwith)
-    coef.re <- coef(xbetw)[tc]
-    vcov.wi <- xwith$vcov
-    vcov.re <- xbetw$vcov[tc,tc]
+      coef.wi <- coef(xwith)
+      coef.re <- coef(xbetw)[tc]
+      vcov.wi <- xwith$vcov
+      vcov.re <- xbetw$vcov[tc,tc]
+
+    } else {
+
+        ## then they are ML
+
+        ## determine which is FE
+        if(is.null(dimnames(x$vcov))) {
+            xwith <- x
+            xbetw <- x2
+        } else {
+            xwith <- x2
+            xbetw <- x
+        }
+            
+
+        tc <- intersect(names(coef(xwith)), names(coef(xbetw)))
+        ## fix because vcov for FE is not named. Aaaargh!
+        wtc <- match(tc, names(coef(xwith)))
+        
+        coef.wi <- coef(xwith)[wtc]
+        coef.re <- coef(xbetw)[tc]
+        vcov.wi <- xwith$vcov[wtc,wtc]
+        vcov.re <- xbetw$vcov[tc,tc]
+        
+    }
+
+
     
     dbeta <- coef.wi - coef.re
     df <- length(dbeta)