Question about estspec.jl


#1

Hi All,

First of all, thank you for all the effort in putting together all this material for Julia, this is extremely useful.

I am experiencing something rather weird. I am using the routines in estspec.jl, trying to replicate the example on the website. I have the following code

T=150
rho=-0.9
e=randn(T)
x=[e[1]]
tmp=e[1]
for t=2:T
    tmp=rho*tmp+e[t]
    push!(x,tmp)
end

w,fw=ar_periodogram(x,"hamming",65)
semilogy(w,fw)

and get the following message

LoadError: DimensionMismatch("arrays could not be broadcast to a common size")

The problem comes from the following line is estspec.jl

I_w = I_w ./ abs(1 - phi .* exp(im.*w)).^2

In fact I_w has dimension (73,) while w has dimension (75,). The problem originates from the convolution in the smooth function.
Surprisingly enough if instead of using a bandwidth of 65, I use 63 or 67, everything works fine, but 61 or 69 fail !
Also note that this occurs with all types of window.

Is there anything wrong I am doing, or something trivial I am missing?

Thank you in advance

Fabrice


#2

Hi Fabrice,

Good to hear from you and sorry for the slow reply.

Thanks for your feedback. It does sound strange. I can’t see anything obvious.

I’ll be able to look into this in more detail next week.

Regards,

John.


#3

Sorry this has sat for so long. I was able to figure it out.

What we are seeing here is another case of “floating point arithmetic is inexact”

The issue can be seen here:

julia> 65/2, round(Int, 65/2)
(32.5,32)

julia> 67/2, round(Int, 67/2)
(33.5,34)

On this line we call round but it just so happens that the computer thinks 65/2 is just slightly less than 32.5 so it rounds down to 32. On the other hand it believes that 67/2 is at least 33.5 so it rounds up to 34.

Changing round to ceil there fixes the problem.

I have made this change and added a test based on the example (see this commit)

Thank you @fabrice.collard for brining this up!


#4

Dear Spencer,

Thanks a lot for your detailed and insightful answer, and rhanks a lot for all your efforts for the community.

Best


#5

Nice detective work.