Em 24 de abril, Allin escreveu:

On Mon, 23 Apr 2012, Henrique Andrade wrote:

I'm trying to make my first Gretl function and, incredibly [:)], I am
getting no success in this endeavor (and I know this is my fault).

The first problem I found is that I can't save the models, so the command
""ARIMA($P,1,$Q)" <- arima $P 1 $Q ; Y --nc" doesn't work and I need to
replace it with "arima $P 1 $Q ; Y --nc".

True, you can't save models by name inside a function. This is part of the "encapsulation" idea for gretl functions: their ability to produce "side effects" is strictly circumscribed. Basically, functions are allowed to print stuff and/or return a value, and nothing else. But if you want to return multiple objects (in a user-friendly form) you can use the bundle data-type for that purpose.

(See also http://www.wfu.edu/~cottrell/tmp/extending.pdf )

This is really great!  But for the time being I will not use the "bundle" because at first I need a better understanding about the construction of simpler functions. But this is now on my "to do list" :-)

Another problem: I can't save the forecasts. The command I'm using is
"fcast Y_hat_$P_1_$Q".

I'm not seeing a problem with that here. This works OK:

function void testarima (series y)
 loop p=1..2
   loop q=1..2
     arima p 1 q ; y --nc
     fcast Y_hat_$p_1_$q
     print Y_hat_$p_1_$q
end function

open fedstl.bin
data exbzus
dataset addobs 12

Of course the series generated by "fcast" are local to the function -- they won't be available outside the function unless you return then in some way.

 So I need to make some modifications in my Hansl code. Please look at this small code:

function list testarima (series y, int p[1:12:2], int q[1:12:2])
    list Forecasts = null
    loop P=1..p
        loop Q=1..q
            arima P 1 Q ; y --nc --quiet
            series dummy_mais = misszero(($uhat >= +2*sd($uhat)))
            series dummy_menos = misszero(($uhat <= -2*sd($uhat)))
            arima P 1 Q ; y dummy_mais dummy_menos --nc --quiet
            matrix teste_t = $coeff./$stderr
            matrix roots = $["roots"]
            if abs(teste_t)>critical(t, $T, 0.025) && abs(roots[,1])>1
                arima P 1 Q ; y --nc
                fcast Y_hat_$P_1_$Q
                list Forecasts += Y_hat_$P_1_$Q
                print Forecasts --byobs
    return Forecasts   
end function

open fedstl.bin
data exbzus exchus
dataset addobs 12
list Lista = testarima(exchus)
print Lista --byobs

The problem now is that my function only prints the forecasts but doesn't store them appropriately inside the list "Projecoes". You can find attached the new version of the AutoARIMA function, but I can show in advance the main modification I introduced:

fcast Y_hat_$P_1_$Q
list Projecoes += Y_hat_$P_1_$Q
print Projecoes --byobs

What I doing wrong now? :-(

Attached you can find two files: "Função AutoARIMA.inp" (with my function), and "Função AutoARIMA.inp" (where you can find the commands that can reproduce my function).

Here are some comments on your function:

1. Using "set halt_on_error off" is a bad idea (perhaps we should get rid of it): if you do that then any errors will cascade and error messages will likely be hard to understand. Do use, instead, the "catch" modifier, and check the $error variable afterwards. This is important for functions that invoke commands such as arima, where failure of convergence is a live possibility.

2. There's no need to "set echo off" or "set messages off" inside a function; that's the default behavior. But you can set those things on for debugging.

3. It's good gretl programming style to use the string representation of loop indices only where a string is actually needed; otherwise just use the numerical value of the index. So:

fcast Y_hat_$P_1_$Q  # OK, strings needed
arima P 1 Q  # strings not needed, don't use "arima $P 1 $Q"

4. Consider limiting the sum of the AR and MA orders to something sane. Your function as written allows up to arima(12, 1, 12), which would surely be grossly over-parameterized.

5. There seem to be some problems with initialization of your arima models that include dummies for observations where the plain arima residuals are greater than two standard deviations. We use nonlinear least squares for initialization of arima with exogenous regressors: that works quite well in some cases but apparently does not work well here. You might consider using "set initvals" to specify your own initialization. For example, in the arima(1,1,1) case:

diff y
ols d_y d_y(-1) x1 x2
matrix m = {$coeff[1], 0.001, $coeff[2], $coeff[3]}
set initvals m
arima 1 1 1 ; y x1 x2 --nc

Thanks a lot for these advices. I will implement them in my code as soon as I can.

Best regards,
Henrique Andrade