[R] Problems with arima0

Pascal Grandeau pgrandeau at free.fr
Fri Dec 28 13:58:06 CET 2001


I use R 1.3.1 and now R 1.4.0 under Windows 98 and 2000 and I have some
problems with arima0.

I made an example : it is just for test (it is zipped in the attached file
modelfr2.zip. The file is in the french csv format with sep= ";" and
dec=".", so I use read.csv2).

I want to do a regression with ARMA (p=1; q=3) errors of the first column
(sendout) on the other (the second column which is the constant is useless),
so this is what I do with R :

> dd<-read.csv2("c:/modelfr2.csv")
> names(dd)     # Just to verify !
 [1] "sendout"         "cste"            "TempLT17"        "lag2"
 [5] "We"              "hivergaz"        "lag7"            "octobre"
 [9] "vendredi"        "Tmoycarre"       "Tmaxi"           "avril"
[13] "TempGT17"        "vacances"        "Lineartrend"     "TempLT171"
[17] "feriesemaine"    "vacoct"          "MoistransT17"    "jeudi"
[21] "mercredi"        "Tmoycarre1"      "WET17"           "vacdec"
[25] "MoistransTnx"    "septembre"       "MoistransDeltaT" "MoistransDTmax"
[29] "Rechauffement"   "samedi"          "janvier"         "automne"
[33] "periode"         "quadratrend"     "aout"            "tempeff"
[37] "MoistransDTmin"  "MoistransTmin"   "hiver"           "juin"

> ari<-arima0(dd[,1],xreg=as.matrix(dd[,-1:-2]),order=c(1,0,3))

> ari$coef
            ar1             ma1             ma2             ma3
intercept
   2.201316e-01    3.206251e-01    3.500737e-02
 2.473141e-02   -5.018293e+06
       TempLT17            lag2              We        hivergaz
lag7
   1.028222e+06    2.917059e-01   -1.566802e+06    2.566452e+06
9.574060e-02
        octobre        vendredi       Tmoycarre           Tmaxi
avril
  -2.295546e+06   -8.999837e+05
 2.736549e+04   -1.100721e+05   -2.612360e+06
       TempGT17        vacances     Lineartrend       TempLT171
feriesemaine
  -1.010776e+06   -1.454308e+05   -1.081656e+07
 2.542789e+05   -9.502527e+05
         vacoct    MoistransT17           jeudi        mercredi
Tmoycarre1
   4.843674e+05
 3.627698e+05   -5.517601e+05   -4.267915e+05   -6.751992e+03
          WET17          vacdec    MoistransTnx       septembre
MoistransDeltaT
  -3.878577e+04   -5.872619e+05    6.298442e+03   -3.339184e+06
3.031612e+05
 MoistransDTmax   Rechauffement          samedi         janvier
automne
  -1.019058e+05    1.479148e+05   -2.047690e+05   -1.078147e+06
4.288875e+05
        periode     quadratrend            aout         tempeff
MoistransDTmin
  -3.476783e+04    1.109006e+07   -2.105259e+05
 3.437851e+05   -6.411536e+04
  MoistransTmin           hiver            juin
   1.043679e+05   -2.434158e+05    1.215601e+05

> ari$code
[1] 0

> diag(ari$var.coef)
            ar1             ma1             ma2             ma3
intercept
   3.132705e-01    3.117503e-01    8.869646e-02
 6.965824e-03   -1.046849e+06
       TempLT17            lag2              We        hivergaz
lag7
   1.171641e+07    9.928982e-05    1.214256e+07    3.898401e+06
9.109221e-05
        octobre        vendredi       Tmoycarre           Tmaxi
avril
   3.532436e+07    3.095100e+06    1.722773e+05    7.896726e+05
2.289079e+07
       TempGT17        vacances     Lineartrend       TempLT171
feriesemaine
   5.674179e+06   -5.124374e+05    9.341290e+05
 4.346287e+06   -4.689245e+06
         vacoct    MoistransT17           jeudi        mercredi
Tmoycarre1
   5.112234e+05    6.510381e+05   -6.286289e+05    2.762435e+06
2.200297e+05
          WET17          vacdec    MoistransTnx       septembre
MoistransDeltaT
   1.511598e+06    1.423986e+06    1.284414e+05    3.616154e+07
6.749349e+05
 MoistransDTmax   Rechauffement          samedi         janvier
automne
   3.733944e+06    4.233166e+06   -1.137565e+06   -1.455328e+06
4.870733e+06
        periode     quadratrend            aout         tempeff
MoistransDTmin
   1.214202e+05   -6.837660e+05   -1.292829e+06    1.243424e+07
2.729197e+07
  MoistransTmin           hiver            juin
   5.775706e+06    6.265115e+06   -9.041021e+06

> ari<-arima0(dd[,1],xreg=as.matrix(dd[,-1:-2]),order=c(1,0,3),delta=-1)

> ari$coef
            ar1             ma1             ma2             ma3
intercept
   2.209774e-01    3.197254e-01    3.456447e-02
 2.639046e-02   -5.018293e+06
       TempLT17            lag2              We        hivergaz
lag7
   1.028222e+06    2.916751e-01   -1.566802e+06    2.566452e+06
9.576419e-02
        octobre        vendredi       Tmoycarre           Tmaxi
avril
  -2.295546e+06   -8.999837e+05
 2.736549e+04   -1.100721e+05   -2.612360e+06
       TempGT17        vacances     Lineartrend       TempLT171
feriesemaine
  -1.010776e+06   -1.454308e+05   -1.081656e+07
 2.542789e+05   -9.502527e+05
         vacoct    MoistransT17           jeudi        mercredi
Tmoycarre1
   4.843674e+05
 3.627698e+05   -5.517601e+05   -4.267915e+05   -6.751992e+03
          WET17          vacdec    MoistransTnx       septembre
MoistransDeltaT
  -3.878577e+04   -5.872619e+05    6.298442e+03   -3.339184e+06
3.031612e+05
 MoistransDTmax   Rechauffement          samedi         janvier
automne
  -1.019058e+05    1.479148e+05   -2.047690e+05   -1.078147e+06
4.288875e+05
        periode     quadratrend            aout         tempeff
MoistransDTmin
  -3.476783e+04    1.109006e+07   -2.105259e+05
 3.437851e+05   -6.411536e+04
  MoistransTmin           hiver            juin
   1.043679e+05   -2.434158e+05    1.215601e+05

> ari$code
[1] 0

> diag(ari$var.coef)
            ar1             ma1             ma2             ma3
intercept
   2.984138e-01    2.971397e-01    8.488563e-02    6.717204e-03
1.976934e+06
       TempLT17            lag2              We        hivergaz
lag7
   7.650781e+05    9.506931e-05    2.345793e+06    9.386527e+05
8.919501e-05
        octobre        vendredi       Tmoycarre           Tmaxi
avril
  -7.568479e+05    3.492818e+05    1.299257e+05   -2.127121e+05
3.733916e+06
       TempGT17        vacances     Lineartrend       TempLT171
feriesemaine
  -9.460373e+05    1.194900e+06    4.093315e+05    9.296041e+05
2.080155e+06
         vacoct    MoistransT17           jeudi        mercredi
Tmoycarre1
   2.358321e+06    1.698256e+06    9.177695e+05    1.282622e+06
1.686617e+05
          WET17          vacdec    MoistransTnx       septembre
MoistransDeltaT
  -4.289835e+04    2.335913e+06    1.251494e+05
 1.733809e+06   -9.864234e+05
 MoistransDTmax   Rechauffement          samedi         janvier
automne
   1.854713e+06    2.666766e+05    4.113436e+04    4.021349e+05
6.950620e+04
        periode     quadratrend            aout         tempeff
MoistransDTmin
   1.000572e+05    1.266920e+06    1.743747e+04    1.754632e+06
2.770144e+06
  MoistransTmin           hiver            juin
   6.929513e+05    6.794158e+05    2.348674e+05
>


So, I get negatives values for diag(ari$var.coef) but the convergence seems
good : ari$code=0 !

I try the same thing with SAS, with the procedure :

PROC ARIMA DATA=WORK._egtemp_;
	IDENTIFY
		VAR=sendout crosscor=(TempLT17 lag2 We hivergaz lag7 octobre vendredi
Tmoycarre Tmaxi avril TempGT17 vacances Lineartrend TempLT171 feriesemaine
vacoct MoistransT17 jeudi mercredi Tmoycarre1 WET17 vacdec MoistransTnx
septembre MoistransDeltaT MoistransDTmax Rechauffement samedi janvier
automne periode quadratrend aout tempeff MoistransDTmin MoistransTmin hiver
juin)
		;
	ESTIMATE
	    	METHOD=ML
		MAXITER=150
		P=(1)
		Q=(1,2,3)
		input=(TempLT17 lag2 We hivergaz lag7 octobre vendredi Tmoycarre Tmaxi
avril TempGT17 vacances Lineartrend TempLT171 feriesemaine vacoct
MoistransT17 jeudi mercredi Tmoycarre1 WET17 vacdec MoistransTnx septembre
MoistransDeltaT MoistransDTmax Rechauffement samedi janvier automne periode
quadratrend aout tempeff MoistransDTmin MoistransTmin hiver juin)
		;
	FORECAST
		ID=Date
		INTERVAL = DAY
		;
RUN;


and here are the results :


Maximum Likelihood Estimation
Parameter	Estimate	Standard Error	t Value	Approx	Lag	Variable	Shift
							Pr>|t|
MU		5677704.9	4178966.3	1.36	0.1743	0	sendout		0
MA1,1		0.13339	0.04068	3.28	0.0010	1	sendout		0
MA1,2		0.31361	0.04518	6.94	<.0001	2	sendout		0
MA1,3		0.12087	0.03985	3.03	0.0024	3	sendout		0
AR1,1		0.98183	0.0090420	108.59	<.0001	1	sendout		0
Intercept	914801.6	31603.0	28.95	<.0001	0	TempLT17	0
NUM2		0.08117	0.02464	3.29	0.0010	0	lag2		0
NUM3		-1260988.9	82649.5	-15.26	<.0001	0	We		0
NUM4		1883728.1	444508.7	4.24	<.0001	0	hivergaz	0
NUM5		0.04916	0.01524	3.23	0.0013	0	lag7		0
NUM6		-1291028.7	678103.9	-1.90	0.0569	0	octobre		0
NUM7		-537376.7	86014.1	-6.25	<.0001	0	vendredi	0
NUM8		20164.3	1676.0	12.03	<.0001	0	Tmoycarre	0
NUM9		-72794.9	14816.2	-4.91	<.0001	0	Tmaxi		0
NUM10		-1735412.6	550555.1	-3.15	0.0016	0	avril		0
NUM11		-777787.6	74365.4	-10.46	<.0001	0	TempGT17	0
NUM12		-176589.8	146977.0	-1.20	0.2296	0	vacances	0
NUM13		-21270415	5072492.7	-4.19	<.0001	0	Lineartrend	0
NUM14		144286.3	29856.4	4.83	<.0001	0	TempLT171	0
NUM15		-746670.3	101013.5	-7.39	<.0001	0	feriesemaine	0
NUM16		143734.4	291312.3	0.49	0.6217	0	vacoct		0
NUM17		197695.0	45790.1	4.32	<.0001	0	MoistransT17	0
NUM18		-180919.2	83234.9	-2.17	0.0297	0	jeudi		0
NUM19		-104867.0	67638.4	-1.55	0.1210	0	mercredi	0
NUM20		2414.5	1067.0	2.26	0.0236	0	Tmoycarre1	0
NUM21		-38796.3	8825.9	-4.40	<.0001	0	WET17		0
NUM22		-890757.8	322913.6	-2.76	0.0058	0	vacdec		0
NUM23		2390.4	1859.2	1.29	0.1985	0	MoistransTnx	0
NUM24		-1688542.1	680933.2	-2.48	0.0131	0	septembre	0
NUM25		193932.8	41955.6	4.62	<.0001	0	MoistransDeltaT	0
NUM26		-59055.8	18549.9	-3.18	0.0015	0	MoistransDTmax	0
NUM27		108044.1	43874.8	2.46	0.0138	0	Rechauffement	0
NUM28		-163420.7	54042.0	-3.02	0.0025	0	samedi		0
NUM29		-231656.4	350368.4	-0.66	0.5085	0	janvier		0
NUM30		662353.6	285263.7	2.32	0.0202	0	automne		0
NUM31		-59281.8	14180.1	-4.18	<.0001	0	periode		0
NUM32		21225993	5095024.2	4.17	<.0001	0	quadratrend	0
NUM33		-120509.9	315387.3	-0.38	0.7024	0	aout		0
NUM34		-65665.9	54566.7	-1.20	0.2288	0	tempeff		0
NUM35		-26494.1	17907.6	-1.48	0.1390	0	MoistransDTmin	0
NUM36		45824.8	46569.2	0.98	0.3251	0	MoistransTmin	0
NUM37		480699.1	289296.2	1.66	0.0966	0	hiver		0
NUM38		210529.4	248285.9	0.85	0.3965	0	juin		0

and the results seems good but are very differents.

I try then another tool which is Matrixer and I got almost the same results
than with SAS.

Finally, with S-Plus 4.5 and arima.mle, I got also almost the same results
than with SAS.

Certainly I did not understood the use of arima0 because something seems to
go wrong.

Another problem is if I put transform.pars=1 then R crash under Windows :
error in ts.dll.

Please can you explain me how I can make a regression with ARMA errors with
R and get
the coefficientd and the estimations of standard errors for the
coefficients.

Maybe I can do this with dse but I did not well understood how to define the
problem with dse.

Thank you very much.

Pascal Grandeau
-------------- next part --------------
A non-text attachment was scrubbed...
Name: modelfr2.zip
Type: application/x-compressed
Size: 31602 bytes
Desc: not available
Url : https://stat.ethz.ch/pipermail/r-help/attachments/20011228/137ffbf6/modelfr2.bin


More information about the R-help mailing list