Ideal Band Pass Filter For Stationary/Non-Stationary Series

For posting your own programs to share with others

Moderators: EViews Gareth, EViews Moderator

souliaris
Posts: 7
Joined: Thu Oct 29, 2009 5:17 am

Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby souliaris » Thu Oct 29, 2009 11:43 am

Implements the IdealBandPass filter as described in Corbae and Ouliaris (2006), "Extracting Cycles From Non Stationary Data", Econometric Theory and Practice, pp. 167-177, Cambridge University Press.
Involves detrending the series in the frequency domain. Series can be stationary, or it can have a unit root.
It tends to produce output similar to the B-K filter, though without the end-point issue (estimates end-points directlty).

Subroutine DFT() computes the discrete fourier transform a series.
Subroutine INVDFT() computed the inverse of the discrete fourier transform. DFT() needs to be called first.
Subroutine DFTE() annihilates appropriate frequencies according to user requirements
Subroutine IdealBandPass() implements the routine (at the end of the code section).

The arguments to this routine are:

'x: series name
's: starting period (can be a fraction of pi, or the number of periods)
'e: ending period (can be a fraction of pi, or the number of periods)
'%s_name = name of the filtered series
'%s_date -- start from this date
'%e_date -- end at this date.

Last two arguments are for series that having missing values at the beginning and end of the workfile range.

'Example call: call IdealBandPass(r_apsp,0.25,0.5,"i_bp_oil","1980m05", "2009m05")

in which: r_apsp is the actual name of the series in the active workfile. "i_bp_oil" will be the name of the filtered series in the same workfile.


Sam Ouliaris, October 29, 2009
souliaris@imf.org



Code: Select all


subroutine dft(series x, series y, string %s_date, string %e_date)

   '   Subroutine to compute the DFT of a (possibly complex) variable.
   '                Y is mostly set to 0.
               '                 Output is saved in two temporary series:  dft_r (real_part), dft_m (imaginary part) 
   '     Sam Ouliaris, October 29, 2009, souliaris@imf.org

   series dft_r = 0.0
   series dft_m = 0.0

   !TwoPi = 2 * 3.141592653589793238462

   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   !p = @dtoo(%s_date)-1

   for !i = 0 to  !n-1

           !ae = -!TwoPi*(!i/!n)
         !rp=0.0
         !cp = 0.0

         for !k = 1 to !n

                                 !a = Cos((!k-1)*!ae)
      !b = Sin((!k-1)*!ae)
                                !rp = !rp + ((x(!k+!p)*!a)-(y(!k+!p)*!b))
      !cp = !cp + ((x(!k+!p)*!b)+(y(!k+!p)*!a))

         next

         dft_r(!i+!p+1) = !rp            
         dft_m(!i+!p+1) = !cp            

   next

   series  dft_r = dft_r/!n
   series dft_m = dft_m/!n

endsub

subroutine invdft(string %s_date, string %e_date)

   'Subroutine to compute the inverse of a DFT.
   'Result is saved in two series dft_inv_r  (real part)
   'dft_inv_m (complex part).

   !TwoPi =  2*3.141592653589793238462

   series dft_inv_r=0.0
   series dft_inv_m=0.0

   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   !p = @dtoo(%s_date)-1

   for !i = 0 to  !n-1

           !ae = !TwoPi*(!i/!n)

         !rp=0.0
         !cp = 0.0

         for !k = 1 to !n

                                 !a = Cos((!k-1)*!ae)
      !b = Sin((!k-1)*!ae)
                                !rp = !rp + ((dft_r(!k+!p)*!a)-(dft_m(!k+!p)*!b))
      !cp = !cp + ((dft_r(!k+!p)*!b)+(dft_m(!k+!p)*!a))

         next

         dft_inv_r(!i+!p+1) = !rp            
         dft_inv_m(!i+!p+1) = !cp            

   next

endsub

subroutine dfte(series x, scalar s, scalar e, string %s_name,string %s_date, string %e_date)

   '   Extracts the relevant frequencies between s, e [inclusive], assuming s, e are between 0 and 1.


   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   series ibp_zero=0.0

   call dft(x,ibp_zero,%s_date,%e_date)

   call rm("ibp_zero")

   !p = @dtoo(%s_date)-1

   !k = 0

   if (@mod(!n,2) = 0.0) then
                    !m = !n / 2
    else
       !m = (!n+1)/2
                endif

   for !j = 1 to !m

      !k = !j/!m
      if ((!k < s) or (!k > e)) then
         if !j = 1 then
            dft_r(!p+1) = 0.0
            dft_m(!p+1) = 0.0
         else
            dft_r(!j+!p)=0.0
            dft_r(!n-(!j-2)+!p)=0.0            
            dft_m(!j+!p)=0.0
            dft_m(!n-(!j-2)+!p)=0.0            
         endif
      endif

   next

   if (s > 0.0) then 'This is to ensure that the 0 frequency is always eliminated for small N

      dft_r(!p+1) = 0.0
      dft_m(!p+1) = 0.0         

   endif

   if (@mod(!n,2)=0.0) then
      if (e < 1.0) then
         dft_r(!m+1+!p) = 0.0
         dft_m(!m+1+p!) = 0.0         
      endif
   endif

   call invdft(%s_date,%e_date) 'Inverse; result goes back to dft_inv_r, dft_inv_e

   series {%s_name} = dft_inv_r

   'Clean up...

   call rm("dft_inv_r")
   call rm("dft_inv_m")
   call rm("dft_r")
   call rm("dft_m")

endsub

subroutine rm(string %s_name)

   if (@isObject(%s_name) = 1) then
      delete {%s_name}
   endif

endsub

' Examplel:   call IdealBandPass(r_gdp,6,32,"Ibp_r_gdp", "start_date", "end_date")

   'r_gdp = name of the series
   '6,32: extract cycle between 6 and 32 periods.  Can also be called with fractions" 1/16=2/32, 1/3=2/6
   

subroutine IdealBandPass(series x, scalar s, scalar e, string %s_name,string %s_date, string %e_date)

   'Implements the IdealBandPass filter,  as in Corbae and Ouliaris (2006), "Extracting Cycles From Non  Stationary Data", Econometric Theory and Practice, pp. 167-177, Cambridge University Press.
                'Series can be stationary, or it can have a unit root.

   'Calling arguments:

      'x: series name
      's: starting period (can be a fraction of pi, or the number of periods)
      'e: ending period (can be a fraction of pi, or the number of periods)
      '%s_name = name of the filtered series
      '%s_date -- start from this date
      '%e_date -- end at this date.
   
   'Last two arguments are for series that having missing values at the beginning and end of the workfile range.

   if (s > 1) then  'Allow the user to use periods when calling the routine.   
                   s1 = 2/e
      e   =  2/s
      s   = 2/e
   endif

   if s <=0 or (s>e) then
                   stop
   endif

   if e <=0 or (e>1) then
                   stop
   endif

   smpl %s_date  %e_date
   series dfft_gg = (@trend+1)/@obs(x) 'Scaled time trend.
   call dfte(dfft_gg,s,e,"dfft_st",%s_date, %e_date) 'Compute the DFT of gg, and extract the relative frequencies.
   call dfte(x,s,e,"dfft_sx",%s_date, %e_date) 'Compute the DFT of x, and extract the relative frequencies.
   equation dft_yy
   freeze(mode=overwrite,dfft_ibp_reg) dft_yy.ls dfft_sx dfft_st 'Run the regression in the time domain.
   series {%s_name} = resid 'Residuals contain the ideal band pass filter.
   call rm("dfft_sx")
   call rm("dfft_st")
   call rm("dfft_gg")
   call rm("dfft_ibp_reg")
   call rm("dft_yy")

endsub

terrya
Posts: 81
Joined: Wed Aug 26, 2009 2:37 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby terrya » Sat Jan 16, 2010 5:05 pm

souliaris wrote:
Implements the IdealBandPass filter as described in Corbae and Ouliaris (2006), "Extracting Cycles From Non Stationary Data", Econometric Theory and Practice, pp. 167-177, Cambridge University Press.
Involves detrending the series in the frequency domain. Series can be stationary, or it can have a unit root.
It tends to produce output similar to the B-K filter, though without the end-point issue (estimates end-points directlty).

Subroutine DFT() computes the discrete fourier transform a series.
Subroutine INVDFT() computed the inverse of the discrete fourier transform. DFT() needs to be called first.
Subroutine DFTE() annihilates appropriate frequencies according to user requirements
Subroutine IdealBandPass() implements the routine (at the end of the code section).

The arguments to this routine are:

'x: series name
's: starting period (can be a fraction of pi, or the number of periods)
'e: ending period (can be a fraction of pi, or the number of periods)
'%s_name = name of the filtered series
'%s_date -- start from this date
'%e_date -- end at this date.

Last two arguments are for series that having missing values at the beginning and end of the workfile range.

'Example call: call IdealBandPass(r_apsp,0.25,0.5,"i_bp_oil","1980m05", "2009m05")

in which: r_apsp is the actual name of the series in the active workfile. "i_bp_oil" will be the name of the filtered series in the same workfile.


Sam Ouliaris, October 29, 2009
souliaris@imf.org



Code: Select all


subroutine dft(series x, series y, string %s_date, string %e_date)

   '   Subroutine to compute the DFT of a (possibly complex) variable.
   '                Y is mostly set to 0.
               '                 Output is saved in two temporary series:  dft_r (real_part), dft_m (imaginary part) 
   '     Sam Ouliaris, October 29, 2009, souliaris@imf.org

   series dft_r = 0.0
   series dft_m = 0.0

   !TwoPi = 2 * 3.141592653589793238462

   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   !p = @dtoo(%s_date)-1

   for !i = 0 to  !n-1

           !ae = -!TwoPi*(!i/!n)
         !rp=0.0
         !cp = 0.0

         for !k = 1 to !n

                                 !a = Cos((!k-1)*!ae)
      !b = Sin((!k-1)*!ae)
                                !rp = !rp + ((x(!k+!p)*!a)-(y(!k+!p)*!b))
      !cp = !cp + ((x(!k+!p)*!b)+(y(!k+!p)*!a))

         next

         dft_r(!i+!p+1) = !rp            
         dft_m(!i+!p+1) = !cp            

   next

   series  dft_r = dft_r/!n
   series dft_m = dft_m/!n

endsub

subroutine invdft(string %s_date, string %e_date)

   'Subroutine to compute the inverse of a DFT.
   'Result is saved in two series dft_inv_r  (real part)
   'dft_inv_m (complex part).

   !TwoPi =  2*3.141592653589793238462

   series dft_inv_r=0.0
   series dft_inv_m=0.0

   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   !p = @dtoo(%s_date)-1

   for !i = 0 to  !n-1

           !ae = !TwoPi*(!i/!n)

         !rp=0.0
         !cp = 0.0

         for !k = 1 to !n

                                 !a = Cos((!k-1)*!ae)
      !b = Sin((!k-1)*!ae)
                                !rp = !rp + ((dft_r(!k+!p)*!a)-(dft_m(!k+!p)*!b))
      !cp = !cp + ((dft_r(!k+!p)*!b)+(dft_m(!k+!p)*!a))

         next

         dft_inv_r(!i+!p+1) = !rp            
         dft_inv_m(!i+!p+1) = !cp            

   next

endsub

subroutine dfte(series x, scalar s, scalar e, string %s_name,string %s_date, string %e_date)

   '   Extracts the relevant frequencies between s, e [inclusive], assuming s, e are between 0 and 1.


   !n = @dtoo(%e_date)-@dtoo(%s_date)+1

   series ibp_zero=0.0

   call dft(x,ibp_zero,%s_date,%e_date)

   call rm("ibp_zero")

   !p = @dtoo(%s_date)-1

   !k = 0

   if (@mod(!n,2) = 0.0) then
                    !m = !n / 2
    else
       !m = (!n+1)/2
                endif

   for !j = 1 to !m

      !k = !j/!m
      if ((!k < s) or (!k > e)) then
         if !j = 1 then
            dft_r(!p+1) = 0.0
            dft_m(!p+1) = 0.0
         else
            dft_r(!j+!p)=0.0
            dft_r(!n-(!j-2)+!p)=0.0            
            dft_m(!j+!p)=0.0
            dft_m(!n-(!j-2)+!p)=0.0            
         endif
      endif

   next

   if (s > 0.0) then 'This is to ensure that the 0 frequency is always eliminated for small N

      dft_r(!p+1) = 0.0
      dft_m(!p+1) = 0.0         

   endif

   if (@mod(!n,2)=0.0) then
      if (e < 1.0) then
         dft_r(!m+1+!p) = 0.0
         dft_m(!m+1+p!) = 0.0         
      endif
   endif

   call invdft(%s_date,%e_date) 'Inverse; result goes back to dft_inv_r, dft_inv_e

   series {%s_name} = dft_inv_r

   'Clean up...

   call rm("dft_inv_r")
   call rm("dft_inv_m")
   call rm("dft_r")
   call rm("dft_m")

endsub

subroutine rm(string %s_name)

   if (@isObject(%s_name) = 1) then
      delete {%s_name}
   endif

endsub

' Examplel:   call IdealBandPass(r_gdp,6,32,"Ibp_r_gdp", "start_date", "end_date")

   'r_gdp = name of the series
   '6,32: extract cycle between 6 and 32 periods.  Can also be called with fractions" 1/16=2/32, 1/3=2/6
   

subroutine IdealBandPass(series x, scalar s, scalar e, string %s_name,string %s_date, string %e_date)

   'Implements the IdealBandPass filter,  as in Corbae and Ouliaris (2006), "Extracting Cycles From Non  Stationary Data", Econometric Theory and Practice, pp. 167-177, Cambridge University Press.
                'Series can be stationary, or it can have a unit root.

   'Calling arguments:

      'x: series name
      's: starting period (can be a fraction of pi, or the number of periods)
      'e: ending period (can be a fraction of pi, or the number of periods)
      '%s_name = name of the filtered series
      '%s_date -- start from this date
      '%e_date -- end at this date.
   
   'Last two arguments are for series that having missing values at the beginning and end of the workfile range.

   if (s > 1) then  'Allow the user to use periods when calling the routine.   
                   s1 = 2/e
      e   =  2/s
      s   = 2/e
   endif

   if s <=0 or (s>e) then
                   stop
   endif

   if e <=0 or (e>1) then
                   stop
   endif

   smpl %s_date  %e_date
   series dfft_gg = (@trend+1)/@obs(x) 'Scaled time trend.
   call dfte(dfft_gg,s,e,"dfft_st",%s_date, %e_date) 'Compute the DFT of gg, and extract the relative frequencies.
   call dfte(x,s,e,"dfft_sx",%s_date, %e_date) 'Compute the DFT of x, and extract the relative frequencies.
   equation dft_yy
   freeze(mode=overwrite,dfft_ibp_reg) dft_yy.ls dfft_sx dfft_st 'Run the regression in the time domain.
   series {%s_name} = resid 'Residuals contain the ideal band pass filter.
   call rm("dfft_sx")
   call rm("dfft_st")
   call rm("dfft_gg")
   call rm("dfft_ibp_reg")
   call rm("dft_yy")

endsub


I just tried this and found that it didn't save the filtered series. I'm using EV 7 but that shouldn't make a difference I think.

EViews Gareth
Fe ddaethom, fe welon, fe amcangyfrifon
Posts: 11252
Joined: Tue Sep 16, 2008 5:38 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby EViews Gareth » Mon Jan 18, 2010 2:19 pm

I think there are a few bugs in the implementation of this code when entering periods, rather than fractions of pi. The code:

Code: Select all


   if (s > 1) then  'Allow the user to use periods when calling the routine.   
                   s1 = 2/e
      e   =  2/s
      s   = 2/e
   endif

is definitely wrong, since s1 has not been declared anywhere. It isn't clear to me what this code is meant to be doing, so I can't quite see how to fix it. I think perhaps the (s>1) condition should also be s>=1, since it seems to me that you should be able to enter s=1 as a valid starting period.


Entering fractions of pi rather than periods should mean that this part of the code is skipped. However the couple of tries I had at entering fractions caused error messages elsewhere to appear...
Follow us on Twitter @IHSEViews

terrya
Posts: 81
Joined: Wed Aug 26, 2009 2:37 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby terrya » Mon Jan 18, 2010 2:57 pm

QMS Gareth wrote:I think there are a few bugs in the implementation of this code when entering periods, rather than fractions of pi. The code:

Code: Select all


   if (s > 1) then  'Allow the user to use periods when calling the routine.   
                   s1 = 2/e
      e   =  2/s
      s   = 2/e
   endif

is definitely wrong, since s1 has not been declared anywhere. It isn't clear to me what this code is meant to be doing, so I can't quite see how to fix it. I think perhaps the (s>1) condition should also be s>=1, since it seems to me that you should be able to enter s=1 as a valid starting period.


Entering fractions of pi rather than periods should mean that this part of the code is skipped. However the couple of tries I had at entering fractions caused error messages elsewhere to appear...



Thanks

souliaris
Posts: 7
Joined: Thu Oct 29, 2009 5:17 am

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby souliaris » Thu Jan 21, 2010 3:02 pm

There were a few bugs in the original code that have (hopefully) been fully addressed in the attached version (fdfilter.prg).

First, as Gareth noted, the initial code checking for "periods" rather than fractions was incorrect (s1 was not declared).
Second, "p!" (invalid syntax) was replaced with "!p" in subroutine dfte().

You can test the revised code on the "Demo.wf1" workfile that is sent out with EViews.

The correct output for the "gdpr" series in the "Demo" dataset is provided in gdp.xls file (attached). The test run (see the last few statements of fdfilter.prg) calculates the FD filter for (6,32) periods -- the usual setting for US real GDP.

Postscript: Gareth pointed out that the s>=1 condition didn't make sense given the subsequent logic for calculating the interval over (0,pi). The correct condition is s>=2. Code has been revised to show s>=2.
Attachments
fdfilter.prg
Eviews code for Ideal Band Pass Filter (Frequency Domain).
(4.56 KiB) Downloaded 911 times
gdp.xls
Expected output from the FD filter for "GDPR", which appears in the "DEMO.WF1" demonstration file sent out with Eviews (Chapter 2).
(33 KiB) Downloaded 547 times
Last edited by souliaris on Thu Jan 21, 2010 3:48 pm, edited 1 time in total.

EViews Gareth
Fe ddaethom, fe welon, fe amcangyfrifon
Posts: 11252
Joined: Tue Sep 16, 2008 5:38 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby EViews Gareth » Thu Jan 21, 2010 3:10 pm

This does seem to be an improvement. I do have one further note though.

I'll admit this is not my area of expertise, so I might be asking silly questions here, but can the start period be equal to 1? I note that you changed your program to allow s>=1, which might indicate that s=1 is a valid input. BUT, if s=1, then e=2/1=2, and you don't allow e>1, so the program still fails...
Follow us on Twitter @IHSEViews

souliaris
Posts: 7
Joined: Thu Oct 29, 2009 5:17 am

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby souliaris » Thu Jan 21, 2010 3:54 pm

Thanks Gareth. Code has been changed in the original post. The correct condition is (s>=2), allowing a maximum e (in terms of fraction of [0 to pi]) of 1. Apologies. Sam

EViews Gareth
Fe ddaethom, fe welon, fe amcangyfrifon
Posts: 11252
Joined: Tue Sep 16, 2008 5:38 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby EViews Gareth » Thu Jan 21, 2010 4:02 pm

Great. For EViews 7 users, I've added a graphical user interface to Sam's program. You'll have to download both my interface program, and Sam's program, since mine is an add-on to his, not a replacement. Note you'll have to remove the two "call" lines at the bottom of Sam's program when running my add-on program.
Attachments
fdfilterui.prg
User Interface for FDFILTER.prg
(2.47 KiB) Downloaded 586 times
Follow us on Twitter @IHSEViews

terrya
Posts: 81
Joined: Wed Aug 26, 2009 2:37 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby terrya » Fri Jan 22, 2010 12:29 am

QMS Gareth wrote:Great. For EViews 7 users, I've added a graphical user interface to Sam's program. You'll have to download both my interface program, and Sam's program, since mine is an add-on to his, not a replacement. Note you'll have to remove the two "call" lines at the bottom of Sam's program when running my add-on program.


I've run the programs and they're great. What would be greater if you could add the GUI to allow a number of series to be filtered in a way that's similar to the Stock prog. I'm doing a study of the results of a number of different filters that involves at least 6 different data series.

EViews Gareth
Fe ddaethom, fe welon, fe amcangyfrifon
Posts: 11252
Joined: Tue Sep 16, 2008 5:38 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby EViews Gareth » Fri Jan 22, 2010 3:27 pm

Simple enough.
Attachments
fdfilterui.prg
(2.69 KiB) Downloaded 574 times
Follow us on Twitter @IHSEViews

terrya
Posts: 81
Joined: Wed Aug 26, 2009 2:37 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby terrya » Fri Jan 22, 2010 4:49 pm

QMS Gareth wrote:Simple enough.


Hi
Just run the program. Stops with error msg: "Syntax error in Series {%S_name} = resids. I was using 2 series.

terrya
Posts: 81
Joined: Wed Aug 26, 2009 2:37 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby terrya » Fri Jan 22, 2010 4:57 pm

terrya wrote:
QMS Gareth wrote:Simple enough.


Hi
Just run the program. Stops with error msg: "Syntax error in Series {%S_name} = resids. I was using 2 series.


Oops! My mistake. Didn't read the instructions properly. Works alright.

randomlearning
Posts: 8
Joined: Thu Jun 03, 2010 5:28 am

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby randomlearning » Tue Jun 15, 2010 10:43 am

Is there a way of using this filter in a loop?

I am trying:

Code: Select all

for %i a b c
call IdealBandPass({%i},0.0625,0.33333,"ibf_{%i}","1950q1","2070q1")
next

but it gives me an error that output name that was specified is reserved.

EViews Gareth
Fe ddaethom, fe welon, fe amcangyfrifon
Posts: 11252
Joined: Tue Sep 16, 2008 5:38 pm

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby EViews Gareth » Tue Jun 15, 2010 10:50 am

try:

Code: Select all

for %i a b c
%temp = "ibf_" + %i
call IdealBandPass({%i},0.0625,0.33333,%temp,"1950q1","2070q1")
next
Follow us on Twitter @IHSEViews

randomlearning
Posts: 8
Joined: Thu Jun 03, 2010 5:28 am

Re: Ideal Band Pass Filter For Stationary/Non-Stationary Series

Postby randomlearning » Tue Jun 15, 2010 11:12 am

Thanks. Now I'm getting the error "Cannot assign string expression to numeric variable in "SERIES {%S_NAME}=resid".


Return to “Program Repository”

Who is online

Users browsing this forum: No registered users and 1 guest