Decouple test dependence from all model parameters (in the namelist).

[charleskawczynski]

Using sampled random numbers for the model input through the namelist means that adding (before others) more sampled random numbers will change the output results, which is confusing. Here's an example:

julia> using Random; Random.seed!(1234)
TaskLocalRNG()

julia> rand(3) # set to namelist["modelA_params"]
3-element Vector{Float64}:
 0.32597672886359486
 0.5490511363155669
 0.21858665481883066

julia> ones(3) # set to namelist["modelB_params"]
3-element Vector{Float64}:
 1.0
 1.0
 1.0

julia> rand(3) # set to namelist["modelC_params"]
3-element Vector{Float64}:
 0.8942454282009883
 0.35311164439921205
 0.39425536741585077

New session:

julia> using Random; Random.seed!(1234)
TaskLocalRNG()

julia> rand(3) # set to namelist["modelA_params"]
3-element Vector{Float64}:
 0.32597672886359486
 0.5490511363155669
 0.21858665481883066

julia> rand(3) # set to namelist["modelB_params"], let's make this random now! Oops, now namelist["modelC_params"] changes
3-element Vector{Float64}:
 0.8942454282009883
 0.35311164439921205
 0.39425536741585077

julia> rand(3) # set to namelist["modelC_params"]
3-element Vector{Float64}:
 0.9531246272848422
 0.7955469475347194
 0.4942498668904206

[mattlevine22]

Do you have a recommendation for a cleaner approach? Currently we set the random seed just before using It to generate the random features, so it is reproducible (unless re-order/insert more random calls after the seed is set). Perhaps we want an accessory function that is imported and called with a random seed and #/style of random features for it to create?

[charleskawczynski]

Do you have a recommendation for a cleaner approach?

Yes, absolutely. Feed in numbers resulting from calibration. This has several advantages:

• Our tests will be more decoupled (the order/insert is exactly the reproducibility issue)
• Our tests should be more robust being that the tuned parameters will be better than using random numbers
• We'll close the loop on actually seeing the profiles for these new models

Is there some reason that we want the model to be able to run (or even test) with random numbers? If you want to add a utility function for reproducibility, a function like

function reproducible_rand(n, set_seed)
    set_seed && Random.seed!(1234)
    return rand(n)
end

would be fine, but we'll still be missing out on the second two points in the above list.

[mattlevine22]

The basic approach of random feature models is to generate a random basis of functions, then fit their coefficients. So, in training, we need to first generate random parameters that form these random bases. We purposely do not calibrate these basis functions. During training, we only coefficients for the random basis functions. Of course, these coefficients are only relevant to the specific random basis realization to which they apply. So, during testing, we want to use the same (randomly generated) bases and their associated learnt coefficients.

On May 26, 2022, at 1:11 PM, Charles Kawczynski @.***> wrote:

Do you have a recommendation for a cleaner approach?

Yes, absolutely. Feed in numbers resulting from calibration. This has several advantages:

Our tests will be more decoupled (the order/insert is exactly the reproducibility issue) Our tests should be more robust being that the tuned parameters will be better than using random numbers We'll close the loop on actually seeing the profiles for these new models Is there some reason that we want the model to be able to run (or even test) with random numbers? If you want to add a utility function for reproducibility, a function like

function reproducible_rand(n, set_seed) set_seed && Random.seed!(1234) return rand(n) end would be fine, but we'll still be missing out on the second two points in the above list.

— Reply to this email directly, view it on GitHub https://github.com/CliMA/TurbulenceConvection.jl/issues/1053#issuecomment-1138920854, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACXVZ355IFZLHYBQQEEPKPDVL7EHNANCNFSM5XB7ATUQ. You are receiving this because you were assigned.

[costachris]

In that case I think we'll need to hard-code the random basis vector to match exactly what was used in the calibration config.

[mattlevine22]

It does get saved in the namelist as rf_fix_ent_params

On May 26, 2022, at 3:25 PM, Costa Christopoulos @.***> wrote:

In that case I think we'll need to hard-code the random basis vector to match exactly what was used in the calibration config.

— Reply to this email directly, view it on GitHub https://github.com/CliMA/TurbulenceConvection.jl/issues/1053#issuecomment-1139074855, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACXVZ362X3JBMDOQOFAP6UDVL7T4TANCNFSM5XB7ATUQ. You are receiving this because you were assigned.

[charleskawczynski]

I see. Then let's add a utility function:

function reproducible_rand(set_seed, seed, ns...)
    set_seed && Random.seed!(seed)
    return rand(ns...)
end

and we can replace all the namelist rand calls with this.

[mattlevine22]

In the Random Features configuration, we make 2 calls to randn and 1 call to rand. When we do Orthogonal Random Features, we build a matrix G which is based on many smaller calls to randn and rand. We can't use the same seed for each call....and it is a bit weird to choose a different seed for every call to rand or randn. Maybe we should build 1 make_RFs utility instead?

[costachris]

Does using reproducible_rand and reproducible_randn in the calibration config in the same way as in (https://github.com/CliMA/TurbulenceConvection.jl/pull/1058 ) solve the problem? After updating the namelist with the calibrated parameters for fno/rfm the only calls to rand/randn will be in rf_fix_ent_params

[charleskawczynski]

Is what's in TC main okay right now? Because we can easily add this (with slight modifications) without affecting CEDMF.

[mattlevine22]

Not sure—what is the question/plan at this point?

Sent from my iPhone

On May 26, 2022, at 5:38 PM, Charles Kawczynski @.***> wrote:

 Is what's in TC main okay right now? Because we can easily add this (with slight modifications) without affecting CEDMF.

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were assigned.

[haakon-e]

Assuming that all random number generation relevant to generation of model parameters is done during namelist construction, would it be an option to generate these numbers using a local seed, not fixing the global seed, e.g.:

# function default_namelist(...)
seed = 1234
rng = MersenneTwister(seed)

namelist["rf_params"] = randn(rng, 5)
...

This does not disturb the global random seed. Consider the following two snippets:

julia> using Random

julia> seed = 1234; rng = MersenneTwister(seed)
MersenneTwister(1234)

julia> Random.seed!(1)
TaskLocalRNG()

julia> rand(3)
3-element Vector{Float64}:
 0.07336635446929285
 0.34924148955718615
 0.6988266836914685

julia> rand(3)
3-element Vector{Float64}:
 0.6282647403425017
 0.9149290036628314
 0.19280811624587546

julia> rand(rng, 3)
3-element Vector{Float64}:
 0.5908446386657102
 0.7667970365022592
 0.5662374165061859

and

julia> using Random

julia> seed = 1234; rng = MersenneTwister(seed)
MersenneTwister(1234)

julia> Random.seed!(1)
TaskLocalRNG()

julia> rand(3)
3-element Vector{Float64}:
 0.07336635446929285
 0.34924148955718615
 0.6988266836914685

julia> rand(rng, 3)
3-element Vector{Float64}:
 0.5908446386657102
 0.7667970365022592
 0.5662374165061859

julia> rand(3)  # numbers generated here identical to 2nd call to rand in the above snippet
3-element Vector{Float64}:
 0.6282647403425017
 0.9149290036628314
 0.19280811624587546

as you can see, the generation of random numbers with the rng seed does not impact the global random seed.

[charleskawczynski]

Not sure—what is the question/plan at this point?

This is what I was thinking: #1062

[mattlevine22]

Ya so haakons use of passing a local RNG for a specific task (eg RF construction) looks like a good way to go.

Basically it just makes the config more robust in cases where we may want to introduce other randomness?

But It doesn’t really change much else right? And currently I think RF is the only time rand and randn are called in the config.

Sent from my iPhone

On May 26, 2022, at 7:03 PM, Charles Kawczynski @.***> wrote:

 Not sure—what is the question/plan at this point?

This is what I was thinking: #1062

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were assigned.

[haakon-e]

I've experienced this "gotcha" a couple of times when running ensembles of my stochastic model. I would strongly prefer if the default behavior of the model is that the "global" seed is not fixed, but that we have an optional knob for testing/reproducibility that fixes the seed for simulations.

Having a local seed for namelist parameters and a global seed for simulations (that is set as needed) would be a good compromise to me.

[costachris]

Isn't passing rng to rand and randn functionally similar to what Charlie proposed with reproducible_rand since a seed is set in that function and calls to it don't impact the global random seed? So adding/removing reproducible_rands doesn't affect the random numbers generated by other calls to the function. Unless I'm missing something

[haakon-e]

I haven't tested Charlie's code, but my understanding is that any call to Random.seed!(seed) sets a global seed, no matter where the code is called.

[charleskawczynski]

In my proposed solution, Random.seed! is never called (if that's what we want in CEDMF).

[mattlevine22]

Haakons is much nicer, because it allows the state of the rng state to evolve. We don’t want to call my_randn(seed=1234) multiple times sequentially, because it will be the same each time.

Another solution is to have an rng be input and output of the rand functions

Sent from my iPhone

On May 27, 2022, at 11:01 AM, Charles Kawczynski @.***> wrote:

 In my proposed solution, Random.seed! is never called (if that's what we want in CEDMF).

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were assigned.

[costachris]

Currently in CEDMF, based on Matt's config, the random feature vector is created like this. What we want to match between TC and CEDMF are the random numbers in rf_fix_ent_params (rf_opt_ent_params are the calibrated parameters that will be hard-coded in the TC namelist after calibration). I think we're just after any function/solution that lets us easily generate the exact same vector of random numbers in the CEMDF config and TC generate_namelist.jl (both of your solutions do). In addition, we don't want the behavior of TC outside of the RFM to change when we do things like change the size of rf_fix_ent_params (which both of your solutions also do AFAIK)

function get_scm_config()
    config = Dict()
    d = 4 # number of Pi groups
    m = 25 # number of random features
    seedval = 2022

    Random.seed!(seedval)
    config["namelist_args"] = [
        ("time_stepping", "dt_min", 1.0),
        ("time_stepping", "dt_max", 2.0),
        ("stats_io", "frequency", 60.0),
        ("grid", "dz", 50.0),
        ("grid", "nz", 80),
        ("turbulence", "EDMF_PrognosticTKE", "area_limiter_power", 0.0),
        ("turbulence", "EDMF_PrognosticTKE", "entr_dim_scale", "inv_z"),
        ("turbulence", "EDMF_PrognosticTKE", "entrainment", "RF"),
        ("turbulence", "EDMF_PrognosticTKE", "entr_pi_subset", ntuple(i -> i, d)),
        (
            "turbulence",
            "EDMF_PrognosticTKE",
            "rf_fix_ent_params",
            vec(cat(2 * pi * rand(2, m, 1), randn(2, m, d), dims = 3)),
        ),
        ("turbulence", "EDMF_PrognosticTKE", "rf_opt_ent_params", vec(cat(randn(2, m), ones(2, d + 1), dims = 2))),
    ]
    return config
end

[haakon-e]

In my proposed solution, Random.seed! is never called (if that's what we want in CEDMF).

Apologies @charleskawczynski, I thought this was your proposed solution?

I see. Then let's add a utility function:

function reproducible_rand(set_seed, seed, ns...)
    set_seed && Random.seed!(seed)
    return rand(ns...)
end

and we can replace all the namelist rand calls with this.

[charleskawczynski]

Apologies @charleskawczynski, I thought this was your proposed solution?

That was the original suggestion. The slightly modified form is in #1062. It adds a reset_seed kwarg (which is instead used in reproducible_rand) to default_namelist, which is set to true. CEDMF can turn it off and it would see the exact same behavior from TC as what's in TC main now.

[charleskawczynski]

I think some of these discussions are getting off-topic from the original intention of the issue, so I'm going to revise the title.

[charleskawczynski]

We can open a new issue about additional features that CEDMF may want, but the solution that I proposed in #1062 is indistinguishable from CEDMF's point of view. Please chime in if there's issues with the solution.

[charleskawczynski]

@haakon-e's suggestion does not actually fix the issue, since the last call to rand(rng, ...) will still depend on how many calls to rand(rng, ...) before it (effectively coupling the results of the tests). Unless we provide a separate rng per call, which @mattlevine22 suggested is awkward--which I agree since it would complicate the interface.

[haakon-e]

How about something like this, @charleskawczynski? edit: updated

# function reproducible_rand(reset_seed, seed, ns...) <- not needed

# ...

function default_namelist(
    case_name::String;
    root::String = ".",
    write::Bool = true,
    set_seed::Bool = false,  # <- I find it much more intuitive that this are false by default
    reset_seed::Bool = false,  # <- ... and this
    seed::Int = 2022,
    truncate_stack_trace::Bool = false,
)

prng = if set_seed
    # I really don't like the idea of setting a global random seed in the namelist, 
    # but can live with it, until I think of something better
    Random.seed!(seed)  
    MersenneTwister(seed)
else
    RandomDevice()
end

param_rng(reset_seed) = reset_seed ? MersenneTwister(seed) : prng

# ...

dims = 5
namelist["rf_params"] = rand(param_rng(reset_seed), dims) 

# ...

end

edit: let me think a little bit more on this ...

[charleskawczynski]

How about something like this

~This seems much more complicated to reason about~. I prefer we enforce reproducibility by default, and allowing CEDMF to opt into stochastic options. The kwargs are there, what's the issue with using them? I'm fine if we want to change to using a rng, but that's besides the point of fixing the issue. So, using something like

function reproducible_rng(set_seed, seed)
    return if set_seed
        Random.seed!(seed)
        MersenneTwister(seed)
    else
        RandomDevice()
    end
end

rng() = reproducible_rng(set_seed, seed)

is fine, but I'd still prefer defaulting set_seed = true.

[mattlevine22]

So I'm a bit lost in this convo---but, I do think Haakon's original suggestion to pass around an RNG STATE is a good one. It allows us to distinguish random components and keep them independently consistent. It would be very easy for me to add this to the config (I believe) without needing accessory functions AND without affecting how anyone else uses randomness.

[charleskawczynski]

Let's discuss this in the next edmf / darpa meeting.

[mattlevine22]

No—my point above is that we may need to do multiple calls to these seeded functions.

Say I need to call these 5-10 times—-then I have to pick a different seed each time. And that is a bit against best use practices with RNGs.

Sent from my iPhone

On May 26, 2022, at 4:30 PM, Costa Christopoulos @.***> wrote:

 Does using reproducible_rand and reproducible_randn in the calibration config in the same way as in (#1058 ) solve the problem?

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were assigned.