Skip to content

gefp() boundaries?

6 messages · Achim Zeileis, bonda

bonda
# 
Hello all,

I have the following two questions: 1) how can I get the values of
boundaries for fluctuation process gefp(), for functionals maxBB, meanL2BB,
etc? 2) how can I get fragments of gefp()-process, e.g., if I have n=200
observations, i=1,2,...,200, and need gefp()[50:100], i.e. from i=50 till
i=100?

Thank you in advance,
Julia

--
View this message in context: http://r.789695.n4.nabble.com/gefp-boundaries-tp3872529p3872529.html
Sent from the R help mailing list archive at Nabble.com.
Achim Zeileis
#
On Tue, 4 Oct 2011, bonda wrote:

            

      
      
      
    

        
Both, maxBB and meanL2BB, are objects of class "efpFunctional" that 
contain several functions including $computeCritval and $boundary that 
compute the critical value c and base boundary b(t), respectively. The 
full boundary is the c * b(t). In case of maxBB and meanL2BB however, the 
boundary is constant anyway, i.e., b(t) = 1, so you essentially need only 
the critical value.

## compute and plot gefp
library("strucchange")
data("durab")
scus <- gefp(y ~ lag, data = durab)
plot(scus, functional = meanL2BB)

## add the critical value (= boundary here) in another color
abline(h = meanL2BB$computeCritval(0.05, nproc = 2), col = 4)

            

      
      
      
    

        
gefp_object$process has the cumulative score process as a "zoo" object 
which preserves the original time scale of the data (if any), e.g.

   proc <- scus$process
   window(proc, start = 2000, end = 2001)
   proc[11:20]

Note, however, that the process has n+1 observations because an additional 
zero is added as the first element (to facilitate visualizations etc.).

See ?gefp, ?efpFunctional for more information and in particular

      Zeileis A. (2006), Implementing a Class of Structural Change
      Tests: An Econometric Computing Approach. _Computational
      Statistics & Data Analysis_, *50*, 2987-3008.
      doi:10.1016/j.csda.2005.07.001.

hth,
Z

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    
1 day later

    

      
      

        


  

        

      bonda
    

    

      
      #
    

  


  

      
Thank you very much for the answer. If I take Poisson model and follow
"Generalized M-fluctuation tests for parameter instability", A. Zeileis and
K. Hornik, Statistica Neerlandica (2007) Vol. 61, N. 4, p. 500-501 (section
4.3):

data("Boston")
n <- 506;
my.X <- as.matrix(cbind(1, Boston["crim"], Boston["age"])); 
my.model <- glm(tax ~ crim + age, family = poisson, data = Boston);
my.psi <- estfun(my.model);
my.mu <- fitted(my.model);
J <- sum(my.mu*my.X%*%t(my.X))/n;
my.process <- apply(as.matrix(my.psi), 2, cumsum)/sqrt(J*n);

gprocess <- gefp(tax ~ crim + age, family = poisson, data=Boston);

then my.process and gprocess$process have to be the same? 

Best regards,
Julia

--
View this message in context: http://r.789695.n4.nabble.com/gefp-boundaries-tp3872529p3878886.html
Sent from the R help mailing list archive at Nabble.com.

  
  


  


      

    

    

      
      

        


  

        

      Achim Zeileis
    

    

      
      #
    

  


  

      
On Thu, 6 Oct 2011, bonda wrote:

        

            

      
      
      
    

        
I guess that this is data("Boston", package = "MASS")?

            

      
      
      
    

        
Applying a Poisson model for count data to a tax rate seems to be very 
awkward.

            

      
      
      
    

        
This should be the covariance matrix, not a scalar.

            

      
      
      
    

        
No. Even if you would have computed the J that you would have been 
interested in (and computed its root correctly), then gefp() still uses a 
different covariance matrix estimator by default.

See ?gefp and the 2006 CSDA paper for more details. Looking at the actual 
code may also help.

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      bonda
    

    

      
      #
    

  


  

      
Well, I am still confused... shouldn't the "made-by-hands"-process (correctly
constructed, of course) and the gefp()-process be similar? For the case of
efp()-process it has been worked, at least... Besides, I'd like to know:
1) if I have the gefp()-process (scaled to unit interval [0;1]), can I get
the unscaled one for the interval [0;n]?
2) is it possible to define a type of the score process like this was in
efp()? I could find that neither in ?gefp, nor in CSDA-paper-2006.
The example was bad and incorrect, the another one:

data("BostonHomicide")

gprocess <- gefp(homicides ~ ahomicides25, family = poisson, data =
BostonHomicide);

my.model <- glm(homicides ~ ahomicides25, family = poisson, data =
BostonHomicide);
J <- vcov(my.model);
J.sqrinv <- solve(sqrtm(J))/n;
my.psi <- estfun(my.model);
my.process <- apply(as.matrix(my.psi), 2, cumsum)%*%t(J.sqrinv)/sqrt(n);

Regards,
Julia





--
View this message in context: http://r.789695.n4.nabble.com/gefp-boundaries-tp3872529p3882197.html
Sent from the R help mailing list archive at Nabble.com.

  
  


  


      

    

    

      
      

        


  

        

      Achim Zeileis
    

    

      
      #
    

  


  

      
On Fri, 7 Oct 2011, bonda wrote:

        

            

      
      
      
    

        
Correctly constructed yes. And if unsure, just read the source code (as 
previously suggested): it's open!

            

      
      
      
    

        
Make your data a "ts" series with the desired index or use the order.by 
argument of gefp().

            

      
      
      
    

        
The 2005 Econ Rev paper on ?gefp explains how various types of processes 
are special cases of the gefp framework by setting suitable functionals.

            

      
      
      
    

        
This is (a) not what the paper says, (b) not on the right scale.

            

      
      
      
    

        
It's really not that difficult. All equations below pertain to the CSDA 
paper:

## fit model
fm <- glm(homicides ~ ahomicides25, family = poisson,
   data = BostonHomicide)
n <- nobs(fm)

## Equation 2
W <- apply(estfun(fm), 2, cumsum) / sqrt(n)

## Equation 3
J <- crossprod(estfun(fm)) / n

## Equation 4
efp1 <- W %*% solve(root.matrix(J))

## via gefp()
efp2 <- gefp(fm, fit = NULL)$process

## using information matrix instead of OPG
efp3 <- W %*% root.matrix(vcov(fm) * n)

## visualize everything
plot(merge(efp2[-1,], efp1, efp3),
   screen = rep(1:2, 3), col = rep(1:3, each = 2))

The green lines pertain to the model with the alternative estimator of the 
covariances and is hence somewhat different.
Z