close Warning: Can't synchronize with repository "(default)" (/var/svn/mms does not appear to be a Subversion repository.). Look in the Trac log for more information.

Changes between Initial Version and Version 1 of Ticket #631


Ignore:
Timestamp:
May 18, 2011, 11:59:22 AM (14 years ago)
Author:
Pedro Gea
Comment:

Legend:

Unmodified
Added
Removed
Modified

  • Ticket #631

    • Property Status changed from new to accepted

  • Ticket #631 – Description

    initial v1  
    44Esto puede ocurrir cuando se ejecuta la previsión de un modelo después de actualizar los datos y no se ha estimado (normalmente, por falta de tiempo). Es necesario adaptar el mecanismo de previsión para tener en cuenta esto y solucionarlo para poder ejecutar la previsión sin tener en cuenta esos datos desconocidos.
    55
    6 Una posible solución puede verse el código que sigue a continuación rescatado de SADD. La función Forecast ejecuta una previsión y utiliza las funciones ActRes y OmitFilterTOL.
    7 {{{
    8 
    9 //////////////////////////////////////////////////////////////////////////////
    10    Serie InterRecur(Set impulsos, Serie resi, Polyn xpiw)
    11 // PURPOSE: Reduce a cero los residuos correspondientes a las fechas de
    12 //          <impulsos> trasmitiendo el valor del residuo segun la funcion de
    13 //          respuesta <xpiw>
    14 //////////////////////////////////////////////////////////////////////////////
    15 {
    16   Date iesima = impulsos[1];
    17   Real valor  = SerDat(resi, iesima);
    18   Serie pulEx = xpiw:(valor * Pulse(iesima, Dating(resi)));
    19   Serie res1  = resi - pulEx;
    20   Set imp1    = impulsos - SetOfDate(iesima);
    21   If(GT(Card(imp1),0), InterRecur(imp1,res1,xpiw), res1)
    22 };
    23 
    24 //////////////////////////////////////////////////////////////////////////////
    25 Serie OmitFilterTOL(Ration piW, Serie res, Serie ind)
    26 //////////////////////////////////////////////////////////////////////////////
    27 {
    28   Polyn xpiW     = Expand(piW,CountS(ind));
    29   Set pulsosF    = Dates(SerTms(ind),First(ind),Last(ind));
    30   Real cardpul   = Card(pulsosF);
    31   If(cardpul, InterRecur(pulsosF,res,xpiW), res)
    32 };
    33 //////////////////////////////////////////////////////////////////////////////
    34 PutDescription(
    35 "funcion que devuelve la serie de correccion de omitidos en la serie de
    36 residuos actualizados.",
    37 OmitFilterTOL);
    38 //////////////////////////////////////////////////////////////////////////////
    39 
    40 //////////////////////////////////////////////////////////////////////////////
    41 Serie ActRes(Serie serRuido, Serie serRes, Real p, Real q, Ration piw)
    42 //////////////////////////////////////////////////////////////////////////////
    43 {
    44   Date  lastRui  = Last(serRuido);
    45   Date  lastRes  = Last(serRes);
    46   Real  lag      = DateDif(Dating(serRuido), lastRes, lastRui);
    47   Real  m        = p + lag;
    48   Serie aIni     = If(q,LastM(serRes,q),CalInd(W,Dating(serRuido)));
    49   Serie noiseIni = LastM(serRuido,m);
    50   If(q,DifEq(piw, noiseIni, aIni), Numerator(piw):noiseIni)
    51 };
    52 //////////////////////////////////////////////////////////////////////////////
    53 PutDescription(
    54 "funcion que actualiza los residuos, a partir de los ordenAR ultimos valores
    55 del ruido y los ordenMA ultimos valores de los residuos. Sus argumentos son:
    56   - serRuido, la serie de ruido del modelo estimado
    57   - serRes, la serie de residuos que se desea actualizar.",
    58 ActRes);
    59 //////////////////////////////////////////////////////////////////////////////
    60 
    61 //////////////////////////////////////////////////////////////////////////////
    62 Set Forecast (Serie original,   // Serie original
    63               Set   boxCox,     // Transformacion Box-Cox
    64               Set   arima,      // Arima Estimado
    65               Serie noise,      // Noise Estimado
    66               Serie filter,     // Filtro Extendido
    67               Serie filterVar,  // SetSum(ForErr(Inputs Estoc)^2)
    68               Serie residuals,  // Residuos Estimados
    69               Real  numFor,     // Numero fde previsiones
    70               Date  origen,     // Origen de prevision (ultimo dato conocido)
    71               Real  alfa,       // 2*(1-DistNormal(1))
    72               Serie cierres)    // Serie de cierres de la serie si los tiene
    73 //////////////////////////////////////////////////////////////////////////////
    74 {
    75   Text  dtnTxt   = "TimeSet "+ DatingName(original) +"; ";
    76 
    77   Polyn   ARI    = ARIMAGetARI(arima);
    78                                  //arima[2]*arima[4];
    79   Polyn   MA     = ARIMAGetMA (arima);//arima[3];//
    80   Real    degARI = Degree(ARI);
    81   Real    degMA  = Degree(MA);
    82   Ration  piw    = ARI/MA;
    83   Ration  psiw   = MA/ARI;
    84 
    85   Date    fsFor  = Succ(origen, Eval(dtnTxt), 1);      // primera prevision
    86   Date    lsFor  = Succ(origen, Eval(dtnTxt), numFor); // ultima  prevision
    87   Real    isEst  = Last(noise) == origen;  // El origen de prevision coincide
    88                                           // con el ultimo dato de estimacion
    89   ////////////////////////////////////////////////////////////////////////////
    90   // Output  y Noise extendido y Residuos actualizados
    91   ////////////////////////////////////////////////////////////////////////////
    92   Serie outputOmit = BoxCoxTransf(original, boxCox);
    93   Serie indOmi     =  // omitidos no estimados
    94     If(isEst, CalInd(W,Eval(dtnTxt)),
    95               SubSer(IsUnknown(outputOmit),
    96                       Succ(Last(noise),Eval(dtnTxt),1),
    97                       origen));
    98 
    99   Serie output   = IfSer(IsUnknown(outputOmit), 0, outputOmit);
    100   Serie filter2For = SubSer(filter, fsFor, lsFor);
    101   Serie noiseExt = noise >> (output - filter);
    102 
    103   Date  lastAvailable    = Min(Last(noiseExt), origen); // ultimo noise
    104   Date  lastAvailableRes = Min(Last(residuals),origen); // ultimo residuo
    105 
    106   Serie noiseNoFil = SubSer(noiseExt, First(noiseExt), lastAvailable);
    107   Serie resEst     = SubSer(residuals, First(residuals), lastAvailableRes);
    108   Real  doAct      = Or((Last(noise) < lastAvailable),
    109                          Last(residuals) < lastAvailableRes);// Ojo con esto
    110   Serie actRes     = If(Not(doAct), resEst, // actualizacion de residuos
    111                             ActRes(noiseNoFil, resEst, degARI, degMA, piw));
    112 
    113   ////////////////////////////////////////////////////////////////////////////
    114   // Residuos y Ruido Filtrado
    115   ////////////////////////////////////////////////////////////////////////////
    116   Date  lsDat     = Last(noiseNoFil);
    117   Serie filtRes   = If(Not(doAct), actRes, OmitFilterTOL(piw, actRes, indOmi));
    118   Serie resActF   = actRes - filtRes;
    119   Serie resActIni = SubSer(CalInd(W,Eval(dtnTxt)),
    120                           Succ(First(resActF), Eval(dtnTxt), -degMA),
    121                           Succ(First(resActF), Eval(dtnTxt), -1));
    122   Serie filNoi   = DifEq(psiw, resActIni<<resActF);
    123   Serie cierre   = SubSer(CalInd(W,Eval(dtnTxt)),
    124                    Succ(lsDat, Eval(dtnTxt), -degARI), lsDat)
    125                    << filNoi;
    126 
    127   Serie hisRes = If(degMA,
    128                     SubSer(filtRes, Succ(lsDat,Eval(dtnTxt),-degMA), lsDat),
    129                     CalInd(W,Dating(filtRes)));
    130 
    131   Serie hisNoi = SubSer((noiseNoFil - cierre),
    132                   Succ(lsDat,Eval(dtnTxt),-degARI), lsDat);
    133 
    134   ////////////////////////////////////////////////////////////////////////////
    135   // Prevision del noise
    136   ////////////////////////////////////////////////////////////////////////////
    137   Serie zero      = SubSer(CalInd(W, Eval(dtnTxt)), fsFor, lsFor);
    138   Serie allRes    = If(degMA,(hisRes>> zero),zero);
    139   Serie prevNoise = SubSer(DifEq(psiw, allRes, hisNoi), fsFor, lsFor);
    140  
    141   Serie prevTrans = filter2For + prevNoise;
    142 
    143   ////////////////////////////////////////////////////////////////////////////
    144   // Error de Prevision
    145   ////////////////////////////////////////////////////////////////////////////
    146   Real  sigma     = StDsS(hisRes);
    147   Serie forErrOut = ErrorForecast(numFor, sigma, origen, psiw, Eval(dtnTxt));
    148  // Set forErrorSet = SetOfSet(ErrorForecast2(numFor, sigma, origen, psiw, Eval(dtnTxt)));
    149   Serie forErr    = SubSer(SqRt(forErrOut^2 + filterVar), fsFor, lsFor);
    150 
    151   ////////////////////////////////////////////////////////////////////////////
    152   // Bandas
    153   ////////////////////////////////////////////////////////////////////////////
    154   Real  numSig  = DistNormalInv(1-(alfa/2));
    155   Serie tSup    = prevTrans + (numSig*forErr);
    156   Serie tInf    = prevTrans - (numSig*forErr);
    157 
    158   ////////////////////////////////////////////////////////////////////////////
    159   // Salida de la funcion
    160   ////////////////////////////////////////////////////////////////////////////
    161   Serie forecast      = BoxCoxInvTransf(prevTrans, boxCox)*Not(cierres);
    162   Serie viewForecast  = original << forecast; 
    163 //  Serie viewForecast  =
    164 //   SubSer(original<<forecast, Last(original), Last(forecast)); 
    165   Serie upperBand     = BoxCoxInvTransf(tSup, boxCox);
    166   Serie lowerBand     = BoxCoxInvTransf(tInf, boxCox);
    167   Set  previsionResults =
    168   [[ original, viewForecast, upperBand, forecast, lowerBand,
    169      prevNoise, filter2For, prevTrans, forErr]];
    170 
    171   Set realSet   = SetOfReal(numFor);
    172   Set aRIMASet  = SetOfPolyn(ARI, MA);
    173   Set dateSet   = SetOfDate(origen, lsFor, lastAvailable, lastAvailableRes);
    174   Set serieSet  = SetOfSerie
    175   (
    176     outputOmit, indOmi, output, noiseExt, noiseNoFil, resEst,
    177     actRes, zero, filtRes, resActF, resActIni, filNoi, cierre,
    178     allRes, hisNoi, forErrOut, forErr, tSup,tInf
    179   );
    180   Set diagnosisSet = SetOfSet(aRIMASet, dateSet, serieSet, realSet);// forErrorSet
    181   previsionResults << [[ diagnosisSet ]]
    182 };
    183 }}}
    184 
    185 
    186 
    187 
     6Una posible solución puede verse el código que se adjunta (sadd_forecast.tol) rescatado de SADD. La función Forecast ejecuta una previsión y utiliza las funciones ActRes y OmitFilterTOL.