0

I am new to Julia (version 1.0.2) and currently trying the @reaction_network from the package DiffEqBiological (also current version, I can't find the version number here):

tspan = (0.0, 50.0);
y0 = [100.0 50.0 0.0 0.0] #[substrate enzyme complex product]
S, E = y0[1], y0[2]
for r in 0.1:0.1:1.0
    println("Creating Michaelis-Menten reaction model...")
    r1, r2, r3 = r, r, r;
    michaelismenten = @reaction_network rType begin
        r1, S + E ⟶ C
        r2, C ⟶ S + E
        r3, C ⟶ P + E
    end

    y = ODEProblem(michaelismenten, y0, tspan)
    sol = solve(y, CVODE_BDF(), reltol=1e-8, abstol=1e-8)
end

When I try to compile this, I get the syntax error: unsupported 'const' declaration on local variable around C:\Users\...\.julia\packages\DiffEqBiological\nujlA\src\reaction_network.jl:447. I tried really hard to find the error and maybe my knowledge is to low to actually understand similar questions, as I tried the answers, but still get the same error message. I tried the usual 

r = 0.1
michaelismenten = @reaction_network rType begin
            r, S + E ⟶ C
            r, C ⟶ S + E
            r, C ⟶ P + E
        end
y = ODEProblem(michaelismenten, y0, tspan)
sol = solve(y, CVODE_BDF(), reltol=1e-8, abstol=1e-8)

and this works perfectly fine. As soon as I add the for-loop to vary my reaction rate, I get this error, meaning this code

for i = 1:10
    r = i/10 # here I thought that Julia maybe did not like the 0.1:0.1:1.0 from above
    michaelismenten = @reaction_network rType begin
                    r, S + E ⟶ C
                    r, C ⟶ S + E
                    r, C ⟶ P + E
                end
end

already throws the same error. That's why, I think that I am doing something wrong with the for loop, that I wanted to use to vary the reaction rate r so that I am able to animate my plots via

anim = @animate for r = 0.1:0.1:1.0
     # create michaelismenten and solve in this for-loop before plotting (see above)
     plot(sol)
end
gif(anim, "$S_$E_RRE.gif", fps = 15)

I hope, that someone can please help me find the solution for my problem. And please help me to make my asked question more clear with formatting.

To summarize the comments a bit:

  • declaring r as global does not work and throws the same error
  • wrapping the michaelismenten reaction network in a function with r as input throws the same error of unsupported 'const' declaration on local variable
Lissa
  • 1
  • 2
  • I don't have the immediate answer to your problem, but it would likely help me, and others, to provide more details. I.e., where are `y0` and `tspan` defined, as well as more details on what versions of Julia and the DiffEqBiological package you are running. Does the first example on https://github.com/JuliaDiffEq/DiffEqBiological.jl work? – jvz Dec 17 '18 at 08:23
  • I added your suggested improvements (definitions of y0, tspan, etc) and tested the first example and this works fine. – Lissa Dec 17 '18 at 08:41
  • I guess it is due to the scoping behaviour. What happens if in your third code block, you put `global r` on the line after starting the `for`-loop? (sorry, don't have access to Julia here to try myself). See also the section "Local scope" in the manual: https://docs.julialang.org/en/v1/manual/variables-and-scoping/index.html – jvz Dec 17 '18 at 08:49
  • It still throws the same error message although I declared the variable as global. I looked into line 447 in reaction_network.jl where the error occurs & it seems to be this function (I don't know how broad I have to look when the error says "around reaction_network.jl:447" `#Turns an array of expressions to a expression block with corresponding expressions. function expr_arr_to_block(exprs) block = :(begin end) foreach(expr -> push!(block.args, expr), exprs) return block end` I am sorry for the weird code formatting, as I tried my best for it to look nice, but it didn't work. – Lissa Dec 17 '18 at 08:57
  • What happens if you wrap the ```michaelismenten = @reaction_network rType begin r, S + E ⟶ C r, C ⟶ S + E r, C ⟶ P + E end``` in a function that takes in `r` and returns `michaelismenten`? If that works, you could just call that function. – jvz Dec 17 '18 at 09:02
  • I tried your suggestion to wrap it in a function that takes `r` as input, but when I compile just the function I get the same error of the unsupported 'const' declaration on a local variable. – Lissa Dec 17 '18 at 09:11

1 Answers1

0

I struggled quite a bit, but this code now works (for those who are interested in a solution ;) ). Thanks, zundertj, for your help! You triggered some ideas in my brain, which lead to this solution!

using Plots
using DifferentialEquations
using Distributions
using DiffEqBiological

function solveRRE(michaelismenten, y0, tspan, p)
    # create and solve problem
    println("Creating ODEProblem...")
    y = ODEProblem(michaelismenten, y0, tspan, p)
    println("Solving Reaction Rate Equation...")
    sol = solve(y, CVODE_BDF(), reltol=1e-8, abstol=1e-8)
    return sol
end

tspan = (0.0, 50.0);
y0 = [100.0 50.0 0.0 0.0] #[substrate enzyme complex product]
S, E = y0[1], y0[2]

michaelismenten = @reaction_network rType begin
    r, S + E ⟶ C
    r, C ⟶ S + E
    r, C ⟶ P + E
end r
anim = @animate for i in 0.1:0.1:1.0
    p = i
    sol = solveRRE(michaelismenten, y0, tspan, p)
    plot(sol)
end
gif(anim, "RRE.gif", fps = 15)

Another possible occuring error is that ffmpeg is not installed on Windows/Linux. For further instructions for that, this website helped me lots: https://video.stackexchange.com/questions/20495/how-do-i-set-up-and-use-ffmpeg-in-windows

Can someone please tell me, how to close this topic as solved? Thanks in advance!

Lissa
  • 1
  • 2