The problem is that the multiplication
y_hih_y = y_hat'*iy_var*y_hat
does not produce a scalar but a 1x1 matrix. You can't assign a 1x1 matrix
to an element that expects a scalar. If you know that the result will be a
1x1 matrix you could extract the scalar by:
y_hih_y = (y_hat'*iy_var*y_hat)[1]
or if you instead declare y_hat as a vector with y_hat = zeros(nvar) you
could write
y_hih_y = dot(y_hat, iy_var*y_hat)
On Wednesday, November 4, 2015 at 3:08:28 AM UTC+1, james...@gmail.com
wrote:
>
> Hi Kristoffer:
>
> Thanks for these suggestions, which are progress on the problem.
>
> The current problem is this error statement.
>
> ERROR: LoadError: MethodError: `convert` has no method matching
> convert(::Type{ForwardDiff.HessianNumber{7,Float64,Tuple{Float64,Float64,Float64,Float64,Float64,Float64,Float64}}},
>
> ::Array{ForwardDiff.HessianNumber{7,Float64,Tuple{Float64,Float64,Float64,Float64,Float64,Float64,Float64}},1})
> This may have arisen from a call to the constructor
> ForwardDiff.HessianNumber{7,Float64,Tuple{Float64,Float64,Float64,Float64,Float64,Float64,Float64}}(...),
> since type constructors fall back to convert methods.
> Closest candidates are:
> ForwardDiff.HessianNumber{N,T,C}(::Any, ::Any)
> call{T}(::Type{T}, ::Any)
> convert{T<:Number}(::Type{T<:Number}, ::Char)
> ...
> in PF_RE_AR1_inner_hhh at
> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/PF_RE_AR1_inner_hhh.jl:48
> in PF_RE_AR1_outer_alt_hhh at
> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/PF_RE_AR1_outer_alt_hhh.jl:79
> in _calc_hessian at
> /home/jim_nason/.julia/v0.4/ForwardDiff/src/api/hessian.jl:98
> in hessian at
> /home/jim_nason/.julia/v0.4/ForwardDiff/src/api/hessian.jl:27
> [inlined code] from
> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/test_RE_AR1_alt.jl:188
> in anonymous at no file:183
> in include at ./boot.jl:261
> in include_from_node1 at ./loading.jl:304
> while loading
> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/test_RE_AR1_alt.jl,
>
> in expression starting on line 131
>
> Attached are revised programs/code for test_RE_AR1_alt.jl,
> PF_RE_AR1_outer_hhh.jl, and PF_RE_AR1_inner_hhh.jl.
>
> ForwardDiff.hessian calls to the last to files. These files implement
> your suggestions and produce the above error statement. The source of the
> error statement seems to be in line 48 of PF_RE_AR1_inner_hhh.jl which is
>
> llhdt[m] = -0.5*(lpi2 + log(det(y_var)) + y_hih_y)
>
> lpi2 is a constant Float64 and llhdt and y_var are
> ForwardDiff.HessianNumber( ....
>
> The problem seems to be that y_hih_y is a 1-element
> Array{ForwardDiff.HessianNumber ....
>
> On line 42 of PF_RE_AR1_inner_hhh.jl
>
> y_hih_y = y_hat'*iy_var*y_hat
>
> where y_hat and iy_var are ForwardDiff.HessianNumber( ....
>
> Suggestions on how to convert y_hih_y to a ForwardDiff.HessianNumber( ....
> are welcome.
>
> Jim
>
>
>
> On Tuesday, November 3, 2015 at 9:12:59 AM UTC-5, Kristoffer Carlsson
> wrote:
>>
>> Just as a side note:
>>
>> Instead of
>> randn(mprt,1)
>>
>> you might just want to use
>>
>> randn(mprt)
>>
>> to create a vector instead of a matrix with one column.
>>
>>
>> On Tuesday, November 3, 2015 at 3:00:37 PM UTC+1, Kristoffer Carlsson
>> wrote:
>>>
>>> The problem is the following:
>>>
>>> When you call PF_RE_AR1_outer_alt from hessian it will send in a
>>> ForwardDiff number as the C_parm argument.
>>>
>>> This line:
>>> sige = C_parm[2] # scale volatility on trend inflation SV
>>>
>>> means that sige will be a ForwardDiff number.
>>>
>>> This line
>>> trd_SV_j1 = randn(mprt,1)*sige + log( trd_SV_0 )
>>>
>>> means that trd_SV_j1 is now a vector of ForwardDiff numbers due to
>>> promotion.
>>>
>>> And the last problem is
>>>
>>> t_SV[:,j] = trd_SV_j1
>>>
>>> where you try to assign a vector of ForwardDiff number to an array which
>>> is supposed to only hold floats.
>>>
>>> The solution to this is to change
>>> t_SV = zeros(mprt, obs)
>>>
>>> to
>>>
>>> t_SV = zeros(eltype(C_parm), mprt, obs)
>>>
>>> for each array that will store an intermediate result. This means that
>>> the element type of the array will always be the same as the element type
>>> of the input argument.
>>>
>>> On Tuesday, November 3, 2015 at 2:43:39 PM UTC+1, james...@gmail.com
>>> wrote:
>>>>
>>>> Hi Kristoffer:
>>>>
>>>> Thank you. Attached are
>>>>
>>>> test_RE_AR_alt.jl
>>>>
>>>> which is the main program that calls to
>>>>
>>>> PF_RE_AR1_outer.jl
>>>>
>>>> The goal of these programs is to estimate a multivariate stochastic
>>>> volatility model using a Bayesian particle Markov Chain Monte Carlo
>>>> simulator.
>>>>
>>>> Using Tony's suggestion solved one problem, but the current issue is
>>>> that running test_RE_AR_alt.jl causes Julia to return
>>>>
>>>>
>>>> LoadError: InexactError()
>>>> in unsafe_setindex! at array.jl:318
>>>> in _unsafe_batchsetindex! at multidimensional.jl:329
>>>> in PF_RE_AR1_outer_alt at
>>>> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/PF_RE_AR1_outer_alt.jl:72
>>>> in _calc_hessian at
>>>> /home/jim_nason/.julia/v0.4/ForwardDiff/src/api/hessian.jl:98
>>>> in hessian at
>>>> /home/jim_nason/.julia/v0.4/ForwardDiff/src/api/hessian.jl:27
>>>> [inlined code] from
>>>> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/test_RE_AR1_alt.jl:182
>>>> in anonymous at no file:177
>>>> in include at ./boot.jl:261
>>>> in include_from_node1 at ./loading.jl:304
>>>> while loading
>>>> /home/jim_nason/jmn_work/smith/NS4/jl_code_Summer2015/RE_rho/MH_PF_test/test_RE_AR1_alt.jl,
>>>>
>>>> in expression starting on line 125
>>>>
>>>>
>>>> Line 72 of PF_RE_AR1_outer_alt.jl is
>>>>
>>>> t_SV[:,j] = trd_SV_j1
>>>>
>>>> where t_SV = zeros(mprt, obs) is allocated prior to entering the loop
>>>> found in PF_RE_AR1_outer.jl, where mprt = 10 and obs = 187 (the length of
>>>> the sample SPF_pi) are constants set in test_RE_AR1_alt.jl and trd_SV_j1
>>>> is
>>>> a random walk process that is mprt x 1 and updated in the loop found in
>>>> PF_RE_AR1_outer.jl. Note that mprt = 10 is to get to the
>>>> ForwardDiff.hessian quickly. More typical is mprt = 500 to 20000.
>>>>
>>>> Since these programs are implementing a simulator, are the indexing
>>>> problems generated by calling to ForwardDiff.hessian problematic?
>>>>
>>>> Also, attached are several more programs related to running
>>>> test_RE_AR1_alt.jl.
>>>>
>>>> Jim
>>>>
>>>>
>>>> On Tuesday, November 3, 2015 at 3:32:20 AM UTC-5, Kristoffer Carlsson
>>>> wrote:
>>>>>
>>>>> As a tip, when you ask for help it is much easier if you post actual
>>>>> runnable code. If you provide a file where "julia thefile.jl" throws an
>>>>> InexactError I'm pretty sure it can be figured out quickly what the
>>>>> problem
>>>>> is.
>>>>>
>>>>> On Tuesday, November 3, 2015 at 12:23:36 AM UTC+1, Tony Kelman wrote:
>>>>>>
>>>>>> It might be a different kind of InexactError then. Try adding a @show
>>>>>> on the rhs value, see what type it is. The array might need to be
>>>>>> constructed in a way that it's capable of holding dual number element
>>>>>> types.
>>>>>>
>>>>>>
>>>>>>
>>>>>> On Mon, Nov 2, 2015 at 3:18 PM -0800, <james...@gmail.com> wrote:
>>>>>>
>>>>>> Hi Tony:
>>>>>>
>>>>>> Thanks, but remember you asked.
>>>>>>
>>>>>> Your first request is answered by
>>>>>>
>>>>>> const dsts = 2 # dim( state vector )
>>>>>>
>>>>>> which is set in the main file that calls to PF_RE_AR1_outer_alt.jl
>>>>>>
>>>>>> In PF_RE_AR1_outer_alt.jl, I set
>>>>>>
>>>>>> cs00V = eye(dsts)
>>>>>> cs00V[2,2] = sigu*sigu/(1.0 - rho1*rho1) # this is the line at
>>>>>> which Julia responds with ERROR: LoadError: InexactError()
>>>>>>
>>>>>> Note that sigu and rho1 are Float64 types.
>>>>>>
>>>>>> If you want, I can send you the complete set of Julia files/functions
>>>>>> generating the problem in ForwardDiff.hessian.
>>>>>>
>>>>>> Best,
>>>>>>
>>>>>> Jim
>>>>>>
>>>>>>
>>>>>> On Monday, November 2, 2015 at 10:31:08 AM UTC-5, Tony Kelman wrote:
>>>>>>
>>>>>> Misread the part where you said dsts = 2. Can you post more of the
>>>>>> code, exactly how is dsts getting passed to these lines?
>>>>>>
>>>>>>
>>>>>> On Monday, November 2, 2015 at 7:28:56 AM UTC-8, Tony Kelman wrote:
>>>>>>
>>>>>> By "the destination array" I meant the line
>>>>>>
>>>>>> cs00V = eye(dsts)
>>>>>>
>>>>>> If dsts is an integer element type there, then cs00V will also have
>>>>>> an integer element type. Try eye(size(dsts)...)
>>>>>>
>>>>>>
>>>>>> On Monday, November 2, 2015 at 5:32:54 AM UTC-8, james...@gmail.com
>>>>>> wrote:
>>>>>>
>>>>>> Hi Tony:
>>>>>>
>>>>>> Thanks for the suggestion.
>>>>>>
>>>>>> Unfortunately, the result is the same.
>>>>>>
>>>>>> I initialized/pre-allocated the destination array for the Hessian as
>>>>>>
>>>>>> hessh = zeros(nparm,nparm)
>>>>>>
>>>>>> Julia's reply from the call to
>>>>>>
>>>>>> hessh = ForwardDiff.hessian(PF_RE_AR1_outer_alt, parms)
>>>>>>
>>>>>> began with
>>>>>>
>>>>>>
>>>>>> ForwardDiff.HessianNumber{7,Float64,Tuple{Float64,Float64,Float64,Float64,Float64,Float64,Float64}}
>>>>>>
>>>>>> ForwardDiff.HessianNumber{7,Float64,Tuple{Float64,Float64,Float64,Float64,Float64,Float64,Float64}}
>>>>>>
>>>>>> followed by the ERROR: LoadError: InexactError() statement.
>>>>>>
>>>>>> Jim
>>>>>>
>>>>>>
>>>>>> On Sunday, November 1, 2015 at 11:35:16 PM UTC-5, Tony Kelman wrote:
>>>>>>
>>>>>> You probably need to initialize the destination array with a floating
>>>>>> point, rather than integer, element type.
>>>>>>
>>>>>>
