repsychling / mixedmodelssim.jl Goto Github PK

Hi,

There's a tutorial for power simulations here from @debruine that was presented during the ZIF research and tutorial workshop and again in the 2020 SMLP Summer School. We found the tutorial very useful but noticed that it does not work with the current version of Julia MixedModelsSim.

We spent some time updating the tutorial and adding some additional subjects that were useful to us. We think our new and updated tutorial may be useful to other users and are preparing a pull request.

We suggest adding the new tutorial to this repo (and not the sim-tutorial repo), because users are more likely to find it here and the sim-tutorial seems to be specific to one presentation. Does this sound good to you, or do you suggest a different place to collect tutorials?

We also have one or two questions that relate directly to our code and the format of the tutorial (ipynb or jmd), that maybe we can discuss in the pull request thread.

• @pianobacki and @lschwetlick

ping: @palday

After we've been accepted into the registry and TagBot has done it's thing, I'll update the README:

• change installation instructions
• add Zenodo
• change repo status to "active"

Phillip told me: "theta is lower cholesky factor of the covariance matrix of the random effects, stored in row-major order".

That is very hard to specify / intuit. Could I specify the covariance matrix and then call lower(cholesky) and vectorize it, or something?

A helper function would be great for this :-)

Since some of the functionality we want to add here (see #1) depends on behavior in MixedModels#master, we should consider tracking the manifest.

A few options

1. more documentation on how to do this directly (since theta is just the scaled standard deviations)
2. add in automatic conversion (with error checking) into update!
3. create zerocorr variants of create_re

I think (2) is the best given the current structuring of the package. If we added a create_theta method, it's less ideal because that method either needs the model as an argument or we need the RE covariance matrices to already be zerocorr'd.

Most users will be using DataFrames anyway, and DataFrames now supports in place joins.

cc @dmbates

Having that done there causes a lot of complexity in trying to get the Project.toml file correct. Can we just build the documentation in the docs directory instead?

Is there are limit on n < 128 = 2^7 for nlevstbl()? If so, could this be raised to 2^8? (Needed to reproduce a simuation). (This might be a good exercise for myself, too, I guess.)

julia> nlevstbl(:subj, 126, :trt => ["C", "T"])
(subj = ["S001", "S002", "S003", "S004", "S005", "S006", "S007", "S008", "S009", "S010"  …  "S117", "S118", "S119", "S120", "S121", "S122", "S123", "S124", "S125", "S126"], trt = ["C", "T", "C", "T", "C", "T", "C", "T", "C", "T"  …  "C", "T", "C", "T", "C", "T", "C", "T", "C", "T"])

julia> nlevstbl(:subj, 128, :trt => ["C", "T"])
ERROR: MethodError: Cannot `convert` an object of type 
  PooledArrays.PooledArray{String{},Int8,1,Array{Int8,1}} to an object of type 
  PooledArrays.PooledArray{String{},Int16,1,Array{Int16,1}}
Closest candidates are:
  convert(::Type{T}, ::Factorization) where T<:AbstractArray at /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/factorization.jl:58
  convert(::Type{T}, ::T) where T<:AbstractArray at abstractarray.jl:16
  convert(::Type{T}, ::T) where T at essentials.jl:218
Stacktrace:
 [1] push!(a::Vector{PooledArrays.PooledVector{String, Int16, Vector{Int16}}}, item::PooledArrays.PooledVector{String, Int8, Vector{Int8}})
   @ Base ./array.jl:972
 [2] nlevstbl(nm::Symbol, n::Int64, vars::Pair{Symbol, Vector{String}})
   @ MixedModelsSim ~/.julia/packages/MixedModelsSim/UXBym/src/utilities.jl:110
 [3] top-level scope
   @ REPL[39]:1

I have no idea why. The documenter workflow is copied almost verbatim from MixedModels.jl and the PR previews work fine.

This issue is used to trigger TagBot; feel free to unsubscribe.

If you haven't already, you should update your TagBot.yml to include issue comment triggers.
Please see this post on Discourse for instructions and more details.

If you'd like for me to do this for you, comment TagBot fix on this issue.
I'll open a PR within a few hours, please be patient!

I think Lisa's design of working with a tidy dataframe from the result of a simulation is better than the current method of returning a NamedTuple of NamedTuples so why not do it directly at the simulation level and eliminate the need to use sim_to_df? I can do it if this seems reasonable. I may need @palday to take a look at the threaded version. It will probably need another lock when pushing results on the end of a Tables.RowTable.