Re: [julia-users] style question about Float initialization
Is there an idiomatic Julian way to deal with these cases without a large
pre-allocation or push? I understand why your suggestion here would be
faster than the original. This application aside, I see a lot of push!()
usage in some of the libs (eg saving state in Optim.jl, readdlm in base,
etc).
Linked lists (they only live in DataStructures.jl?) have all the obvious
problems that Stefan points out (scattered in memory), but prevent big
re-allocations if you can take the performance hit on any calculations.
Happy to be schooled here.
Cameron
On Thursday, May 1, 2014, John Myles White wrote:
> For this kind of problem, you could allocate a big matrix in advance and
> fill it in column-by-column, rather than push onto an array. Something
> roughly like the following would work:
>
> theta = Array(Float64, 3, 10_000)
> theta[:, 1] = [2.0, 3.0, 4.0]
> for k= 2:10_000
> theta_cur = theta[:, k - 1]
> theta_prop = theta_cur + rand(3)
> likelihoodratio = computelr(theta_prop, theta_curr)
> if rand() < min(1,likelihoodratio)
> theta[:, k] = theta_prop
> else
> theta[:, k] = theta_curr
> end
> end
>
> You could make something much faster still with devectorization. The code
> for simulated annealing in Optim.jl might be helpful in this regard,
> although it doesn't retain the full history of the chain.
>
> -- John
>
> On May 1, 2014, at 9:46 AM, Ethan Anderes
> >
> wrote:
>
> > Ok, here is quickest example that I came up with.
> >
> >
> > theta_init = [2.0, 3.0, 4.0]
> > chain = Array{Float64,1}[init]
> > for k=1:10_000
> > theta_prop = chain[end]
> > theta_prop += rand(3)
> > likelihoodratio = computelr(theta_prop, theta_curr)
> > if rand() < minimum([1,likelihoodratio])
> > push!(chain, theta_prop)
> > else
> > push!(chain, theta_curr)
> > end
> > end
> >
> > I'm generating a markov chain of vectors. The code proposes a new state
> vector, theta_prop, and pushes it to chain if it satisfies a criterion. The
> problem is that I may run it for 10_000 steps, look at the results and then
> decide to run it for longer... so preinitializing isn't as natural.
> >
> > BTW: thanks a ton...I can't believe I'm having Stefan look at my code!
> >
> >
>
>
Re: [julia-users] style question about Float initialization
It's not terrible, and you can help it out with sizehint, but it's just possible that it will have to move the data around a bunch. > On May 1, 2014, at 1:06 PM, Ethan Anderes wrote: > > Ok, so I'm getting the message that my push! stuff isn't a good default > goto. I just felt so natural but I can see that John's code is much more > performance critical. I'll wean myself off push! :)
Re: [julia-users] style question about Float initialization
2.0 is a Float64 everywhere – it's only integer literals whose type depend on the system. Even in old 32-bit machines, the FPUs had 64-bit floats, so there was never a place where that was an issue; but you really don't want to be using 64-bit integers on a 32-bit machine. > On May 1, 2014, at 1:01 PM, Ethan Anderes wrote: > > Jameson: > > My response is going to expose my ignorance here ... but I wasn't worried > that Float64 wouldn't be fast, but rather that on someones machine > typeof(2.0) would say Float32 instead of Float64. In which case my code > wouldn't run without changing all the Float64 type initialization to Float32. >
Re: [julia-users] style question about Float initialization
Ok, so I'm getting the message that my push! stuff isn't a good default goto. I just felt so natural but I can see that John's code is much more performance critical. I'll wean myself off push! :)
Re: [julia-users] style question about Float initialization
> On May 1, 2014, at 12:46 PM, Ethan Anderes wrote: > > BTW: thanks a ton...I can't believe I'm having Stefan look at my code! Well, I can't believe someone can't believe I'm looking at their code, so that makes us even ;-)
Re: [julia-users] style question about Float initialization
Jameson: My response is going to expose my ignorance here ... but I wasn't worried that Float64 wouldn't be fast, but rather that on someones machine typeof(2.0) would say Float32 instead of Float64. In which case my code wouldn't run without changing all the Float64 type initialization to Float32.
Re: [julia-users] style question about Float initialization
For this kind of problem, you could allocate a big matrix in advance and fill
it in column-by-column, rather than push onto an array. Something roughly like
the following would work:
theta = Array(Float64, 3, 10_000)
theta[:, 1] = [2.0, 3.0, 4.0]
for k= 2:10_000
theta_cur = theta[:, k - 1]
theta_prop = theta_cur + rand(3)
likelihoodratio = computelr(theta_prop, theta_curr)
if rand() < min(1,likelihoodratio)
theta[:, k] = theta_prop
else
theta[:, k] = theta_curr
end
end
You could make something much faster still with devectorization. The code for
simulated annealing in Optim.jl might be helpful in this regard, although it
doesn't retain the full history of the chain.
-- John
On May 1, 2014, at 9:46 AM, Ethan Anderes wrote:
> Ok, here is quickest example that I came up with.
>
>
> theta_init = [2.0, 3.0, 4.0]
> chain = Array{Float64,1}[init]
> for k=1:10_000
> theta_prop = chain[end]
> theta_prop += rand(3)
> likelihoodratio = computelr(theta_prop, theta_curr)
> if rand() < minimum([1,likelihoodratio])
> push!(chain, theta_prop)
> else
> push!(chain, theta_curr)
> end
> end
>
> I'm generating a markov chain of vectors. The code proposes a new state
> vector, theta_prop, and pushes it to chain if it satisfies a criterion. The
> problem is that I may run it for 10_000 steps, look at the results and then
> decide to run it for longer... so preinitializing isn't as natural.
>
> BTW: thanks a ton...I can't believe I'm having Stefan look at my code!
>
>
Re: [julia-users] style question about Float initialization
Cameron: I'm not sure what linked lists in Julia are but I will look into it. Thanks!
Re: [julia-users] style question about Float initialization
I don't know of a machine that can run Julia for which Float64 is not fast.
However, you could make a typealias and use that for embedded processors
without an fpu for the eventual day when someone ports Julia.
On Thursday, May 1, 2014, Ethan Anderes wrote:
> Oop the second line of the code block should read
>
> chain = Array{Float64,1}[theta_init]
>
Re: [julia-users] style question about Float initialization
Stefan:
Ok so maybe that wasn't the best example to illustrate the Float64 issue. I
guess using the parametric types in functions works so solve this (but not the
push! issue)
function chaininit{T}(theta_init::Vector{T}, n)
chain = Array{T,1}[theta_init]
for k=1:n
theta_prop = chain[end]
theta_prop += rand(3)
likelihoodratio = computelr(theta_prop, theta_curr)
if rand() < minimum([1,likelihoodratio])
push!(chain, theta_prop)
else
push!(chain, theta_curr)
end
end
chain
end
function chainadd!{T}(chain::Vector{Vector{T}}, n)
for k=1:n
theta_prop = chain[end]
theta_prop += rand(3)
likelihoodratio = computelr(theta_prop, theta_curr)
if rand() < minimum([1,likelihoodratio])
push!(chain, theta_prop)
else
push!(chain, theta_curr)
end
end
end
Re: [julia-users] style question about Float initialization
Would linked lists be a more natural storage for these events?
On Thursday, May 1, 2014, Ethan Anderes wrote:
> Oop the second line of the code block should read
>
> chain = Array{Float64,1}[theta_init]
>
Re: [julia-users] style question about Float initialization
Oop the second line of the code block should read
chain = Array{Float64,1}[theta_init]
Re: [julia-users] style question about Float initialization
Ok, here is quickest example that I came up with.
theta_init = [2.0, 3.0, 4.0]
chain = Array{Float64,1}[init]
for k=1:10_000
theta_prop = chain[end]
theta_prop += rand(3)
likelihoodratio = computelr(theta_prop, theta_curr)
if rand() < minimum([1,likelihoodratio])
push!(chain, theta_prop)
else
push!(chain, theta_curr)
end
end
I'm generating a markov chain of vectors. The code proposes a new state vector,
theta_prop, and pushes it to chain if it satisfies a criterion. The problem is
that I may run it for 10_000 steps, look at the results and then decide to run
it for longer... so preinitializing isn't as natural.
BTW: thanks a ton...I can't believe I'm having Stefan look at my code!
Re: [julia-users] style question about Float initialization
If it's in a function, you can just make the function generic, but a bigger
issue is that that push! pattern is bound to be quite slow. It's *much*
better if you can allocate the entire array up front. I'm having a hard
time coming up with examples based on what you wrote – it's a little too
vague. Do you have a more specific example?
On Thu, May 1, 2014 at 12:19 PM, Ethan Anderes wrote:
> I have a question which is somewhat related to this question (
> https://groups.google.com/forum/m/?fromgroups#!topic/julia-users/bldp27KJCjY)
> but mine is more of a style question rather than a language design question.
>
> I find myself often filling an array of Array{Float64,1}'s in julia.
> I usually do this by initializing an empty array then using push!
>
> x = Array{Float64,1}[]
> while someconditional
> nextvec = getvec(...)
> push!(x,nextvec)
> end
>
> I'm wondering how to make the Float64 specification more flexible so that
> the code works when used by someone on a machine using Float32, for
> example, but still keep the performance of having the entries of x all a
> homogeneous type. I would be tempted to initialize
> x=Array{FloatingPoint,1}[] but that allows some entries of x to be
> Array{Float64,1} and others to be Array{Float32,1} for example.
>
>
>
>
