reworkhow / JWAS.jl

Julia for Whole-genome Analysis Software
http://QTL.rocks
GNU General Public License v2.0
96 stars 44 forks source link

ERROR: MethodError: no method matching zero(::Type{Any}) when fitting multivariate Bayesian genomic regression #64

Closed markcharder closed 4 years ago

markcharder commented 4 years ago

I am hoping to fit a multiple trait Bayesian genomic prediction model using the JWAS package. I've managed to run the tutorials without any issues. However, when I apply the model to my own data, I get an error message. I've checked and re-checked the format of the data, which, as far as I can tell, is fine. In fact, the analysis will run for a single trait model but only throws the error for a multi-trait model.

Here is a snippet of some preliminary code that I ran recently:

using JWAS, CSV, DataFrames
phenotypes = CSV.read("phenotypes_filtered-jwas.csv", delim=",", header=true)
model_equation = "y3 = intercept + ID
                  y2 = intercept + ID"
model = build_model(model_equation)
add_genotypes(model,"genotypes_filtered-jwas.csv", separator=',', header=true)
out=runMCMC(model,phenotypes,methods="BayesC",output_samples_frequency=100, estimatePi=true)

This is what the phenotypes file looks like:

218×4 DataFrame
│ Row │ ID     │ y1      │ y2      │ y3      │
│     │ String │ Float64 │ Float64 │ Float64 │
├─────┼────────┼─────────┼─────────┼─────────┤
│ 1   │ a53    │ 387.203 │ 11.1667 │ 1127.78 │
│ 2   │ a48    │ 1832.58 │ 12.6667 │ 1814.44 │
│ 3   │ a54    │ 1329.98 │ 12.0    │ 1641.95 │
│ 4   │ a49    │ 1095.57 │ 13.4    │ 1352.56 │
│ 5   │ a57    │ 1896.0  │ 8.83333 │ 2248.22 │
│ 6   │ a58    │ 2461.06 │ 13.6667 │ 2436.69 │
│ 7   │ a52    │ 1852.94 │ 13.6    │ 1834.59 │
│ 8   │ a55    │ 393.055 │ 13.0    │ 1144.82 │
│ 9   │ a50    │ 1328.05 │ 17.5    │ 1314.9  │
│ 10  │ a59    │ 1635.03 │ 14.3333 │ 1618.84 │
│ 11  │ a51    │ 1334.64 │ 18.1667 │ 1321.43 │
│ 12  │ a60    │ 592.382 │ 11.8333 │ 1725.38 │
│ 13  │ a39    │ 1559.66 │ 14.8    │ 1544.21 │
│ 14  │ a44    │ 1277.02 │ 10.3333 │ 1887.23 │
│ 15  │ a29    │ 780.165 │ 16.1667 │ 925.097 │
⋮
│ 203 │ a126   │ 995.957 │ 13.0    │ 1180.98 │
│ 204 │ a211   │ 1052.86 │ 15.0    │ 1299.83 │
│ 205 │ a125   │ 704.754 │ 12.8333 │ 1381.87 │
│ 206 │ a113   │ 2502.71 │ 7.66667 │ 2477.93 │
│ 207 │ a120   │ 1288.7  │ 13.5    │ 1528.1  │
│ 208 │ a228   │ 3603.53 │ 7.83333 │ 3567.85 │
│ 209 │ a135   │ 2205.77 │ 14.8333 │ 2615.53 │
│ 210 │ a129   │ 1293.42 │ 12.5    │ 1911.46 │
│ 211 │ a132   │ 132.957 │ 16.0    │ 633.127 │
│ 212 │ a140   │ 598.38  │ 12.6667 │ 884.306 │
│ 213 │ a155   │ 527.888 │ 17.3333 │ 1035.07 │
│ 214 │ a156   │ 388.525 │ 18.1667 │ 761.813 │
│ 215 │ a177   │ 2365.01 │ 19.8333 │ 2341.59 │
│ 216 │ a190   │ 107.565 │ 17.6667 │ 313.297 │
│ 217 │ a210   │ 1107.33 │ 10.5    │ 1636.44 │
│ 218 │ a246   │ 1785.46 │ 17.2    │ 1767.78 │

This is what the genotypes file looks like:

218×13366 DataFrame. Omitted printing of 13346 columns
│ Row │ ID     │ m1    │ m2    │ m3    │ m4    │ m5    │ m6    │ m7    │ m8    │ m9    │ m10   │ m11   │ m12   │ m13   │ m14   │ m15   │ m16   │ m17   │ m18   │ m19   │
│     │ String │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │ Int64 │
├─────┼────────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┼───────┤
│ 1   │ a53    │ 2     │ 2     │ 1     │ 1     │ 1     │ 1     │ 0     │ 2     │ 0     │ 1     │ 1     │ 1     │ 1     │ 1     │ 1     │ 0     │ 2     │ 2     │ 0     │
│ 2   │ a48    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 3   │ a54    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 4   │ a49    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 5   │ a57    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 6   │ a58    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 7   │ a52    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 8   │ a55    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 9   │ a50    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 10  │ a59    │ 2     │ 2     │ 1     │ 1     │ 1     │ 1     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 11  │ a51    │ 2     │ 2     │ 1     │ 0     │ 0     │ 2     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 12  │ a60    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 13  │ a39    │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 14  │ a44    │ 2     │ 2     │ 1     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 15  │ a29    │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 0     │ 2     │ 2     │ 0     │ 2     │ 0     │ 1     │ 0     │ 0     │ 0     │ 2     │ 2     │ 2     │
⋮
│ 203 │ a126   │ 2     │ 2     │ 1     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 1     │ 0     │ 2     │ 2     │
│ 204 │ a211   │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 0     │
│ 205 │ a125   │ 2     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │
│ 206 │ a113   │ 2     │ 0     │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 2     │ 0     │
│ 207 │ a120   │ 2     │ 2     │ 0     │ 2     │ 0     │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │
│ 208 │ a228   │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 209 │ a135   │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 210 │ a129   │ 1     │ 1     │ 0     │ 1     │ 2     │ 0     │ 1     │ 2     │ 1     │ 1     │ 1     │ 1     │ 2     │ 1     │ 1     │ 1     │ 1     │ 1     │ 1     │
│ 211 │ a132   │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 1     │ 0     │ 0     │ 0     │ 2     │ 2     │ 1     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │
│ 212 │ a140   │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │
│ 213 │ a155   │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 2     │ 0     │ 2     │ 0     │ 0     │ 1     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │
│ 214 │ a156   │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │
│ 215 │ a177   │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │
│ 216 │ a190   │ 2     │ 0     │ 2     │ 2     │ 0     │ 2     │ 1     │ 0     │ 0     │ 0     │ 2     │ 2     │ 1     │ 2     │ 2     │ 2     │ 2     │ 2     │ 0     │
│ 217 │ a210   │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 1     │ 2     │ 2     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 2     │ 0     │ 1     │ 2     │
│ 218 │ a246   │ 2     │ 2     │ 0     │ 2     │ 2     │ 0     │ 0     │ 2     │ 0     │ 0     │ 2     │ 2     │ 2     │ 0     │ 0     │ 0     │ 2     │ 2     │ 0     │

The error message that I get from running this model is as follows:

ERROR: MethodError: no method matching zero(::Type{Any})
Closest candidates are:
  zero(::Type{Union{Missing, T}}) where T at missing.jl:105
  zero(::Type{Missing}) at missing.jl:103
  zero(::Type{LibGit2.GitHash}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/LibGit2/src/oid.jl:220
  ...
Stacktrace:
 [1] zero(::Type{Any}) at ./missing.jl:105
 [2] _spgetindex at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/sparsevector.jl:777 [inlined]
 [3] getindex at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/sparsevector.jl:782 [inlined]
 [4] dot(::SparseArrays.SparseVector{Any,Int64}, ::Array{Float32,1}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/LinearAlgebra/src/generic.jl:896
 [5] * at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/LinearAlgebra/src/adjtrans.jl:245 [inlined]
 [6] Gibbs(::SparseArrays.SparseMatrixCSC{Any,Int64}, ::Array{Float32,1}, ::SparseArrays.SparseMatrixCSC{Any,Int64}) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/buildMME/solver.jl:151
 [7] macro expansion at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/MCMC/MT_MCMC_BayesianAlphabet.jl:143 [inlined]
 [8] macro expansion at /home/mark/.julia/packages/ProgressMeter/g1lse/src/ProgressMeter.jl:717 [inlined]
 [9] MT_MCMC_BayesianAlphabet(::Int64, ::JWAS.MME, ::DataFrame; Rinv::Array{Float32,1}, burnin::Int64, Pi::Dict{Array{Float64,1},Float64}, estimatePi::Bool, estimate_variance::Bool, estimateScale::Bool, sol::Array{Float32,1}, outFreq::Int64, output_samples_frequency::Int64, methods::String, missing_phenotypes::Bool, constraint::Bool, update_priors_frequency::Int64, output_file::String, causal_structure::Bool) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/MCMC/MT_MCMC_BayesianAlphabet.jl:139
 [10] runMCMC(::JWAS.MME, ::DataFrame; heterogeneous_residuals::Bool, chain_length::Int64, starting_value::Bool, burnin::Int64, output_samples_frequency::Int64, update_priors_frequency::Int64, methods::String, estimate_variance::Bool, Pi::Float64, estimatePi::Bool, estimateScale::Bool, single_step_analysis::Bool, pedigree::Bool, categorical_trait::Bool, missing_phenotypes::Bool, constraint::Bool, causal_structure::Bool, outputEBV::Bool, output_heritability::Bool, seed::Bool, printout_model_info::Bool, printout_frequency::Int64, big_memory::Bool, double_precision::Bool, output_samples_file::String, output_samples_for_all_parameters::Bool) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/JWAS.jl:286
 [11] top-level scope at REPL[11]:1

I would really appreciate if anybody can figure out what the issue might be. Sorry for such a long message.

Mark

markcharder commented 4 years ago

Just a quick update to this issue. I've found that when I run the example multivariate model found here), I stumble upon the same issue. For example:

using JWAS,JWAS.Datasets,DataFrames,CSV

phenofile  = Datasets.dataset("example","phenotypes.txt")
pedfile    = Datasets.dataset("example","pedigree.txt")
genofile   = Datasets.dataset("example","genotypes.txt")

phenotypes = CSV.read(phenofile,delim = ',',header=true,missingstrings=["NA"])
pedigree   = get_pedigree(pedfile,separator=",",header=true);

model_equation2 ="y1 = intercept + x1 + x3 + ID + dam
                  y2 = intercept + x1 + x2 + x3 + ID
                  y3 = intercept + x1 + x1*x3 + x2 + ID";

model2 = build_model(model_equation2);

set_covariate(model2,"x1");

set_random(model2,"x2");
set_random(model2,"ID dam",pedigree);

add_genotypes(model2,genofile,separator=',');

out2=runMCMC(model2,phenotypes,methods="BayesC",estimatePi=true);

This produces the same error:


ERROR: MethodError: no method matching zero(::Type{Any})
Closest candidates are:
  zero(::Type{Union{Missing, T}}) where T at missing.jl:105
  zero(::Type{Missing}) at missing.jl:103
  zero(::Type{LibGit2.GitHash}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/LibGit2/src/oid.jl:220
  ...
Stacktrace:
 [1] zero(::Type{Any}) at ./missing.jl:105
 [2] _zeros_eltypes at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/higherorderfns.jl:204 [inlined]
 [3] _noshapecheck_map(::typeof(+), ::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float64,Int64}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/higherorderfns.jl:159
 [4] map(::typeof(+), ::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float64,Int64}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/higherorderfns.jl:1153
 [5] +(::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float64,Int64}) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.4/SparseArrays/src/sparsematrix.jl:1648
 [6] addVinv(::JWAS.MME) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/buildMME/random_effects.jl:209
 [7] getMME(::JWAS.MME, ::DataFrame) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/buildMME/build_MME.jl:261
 [8] init_mixed_model_equations(::JWAS.MME, ::DataFrame, ::Bool) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/JWAS.jl:581
 [9] runMCMC(::JWAS.MME, ::DataFrame; heterogeneous_residuals::Bool, chain_length::Int64, starting_value::Bool, burnin::Int64, output_samples_frequency::Int64, update_priors_frequency::Int64, methods::String, estimate_variance::Bool, Pi::Float64, estimatePi::Bool, estimateScale::Bool, single_step_analysis::Bool, pedigree::Bool, categorical_trait::Bool, missing_phenotypes::Bool, constraint::Bool, causal_structure::Bool, outputEBV::Bool, output_heritability::Bool, seed::Bool, printout_model_info::Bool, printout_frequency::Int64, big_memory::Bool, double_precision::Bool, output_samples_file::String, output_samples_for_all_parameters::Bool) at /home/mark/.julia/packages/JWAS/4xynz/src/1.JWAS/src/JWAS.jl:252
 [10] top-level scope at REPL[19]:1
fmorgante commented 4 years ago

I get the same error with Julia 1.4.1, but it works as expected with Julia 1.1.1.

rohanLuigi commented 4 years ago

I ran this example successfully with Julia 1.3.1

reworkhow commented 4 years ago

@markcharder @fmorgante Thanks for letting me know. There was a bug with Julia 1.4. I have released JWAS 0.8.4. It should work in Julia1.4 now.

More details about the bug can be found here. https://github.com/reworkhow/JWAS.jl/pull/59

fmorgante commented 4 years ago

@reworkhow thank you for fixing this! I also noticed that when setting center=false in add_genotypes() I get the following error (using the multivariate example above):

julia> add_genotypes(model2,genofile,separator=',', center=false);
The delimiter in genotypes.txt is ','.
The header (marker IDs) is provided in genotypes.txt.
ERROR: MethodError: objects of type Bool are not callable
Stacktrace:
 [1] #readgenotypes#52(::Char, ::Bool, ::Bool, ::Function, ::String) at /home/fmorgante/.julia/packages/JWAS/fZ8iA/src/1.JWAS/src/markers/readgenotypes.jl:108
 [2] #readgenotypes at /home/fmorgante/.julia/packages/JWAS/fZ8iA/src/1.JWAS/src/markers/readgenotypes.jl:0 [inlined]
 [3] #add_genotypes#51(::Char, ::Bool, ::Bool, ::Bool, ::Int64, ::Function, ::JWAS.MME, ::String, ::Bool) at /home/fmorgante/.julia/packages/JWAS/fZ8iA/src/1.JWAS/src/markers/readgenotypes.jl:61
 [4] #add_genotypes at ./none:0 [inlined] (repeats 2 times)
 [5] top-level scope at none:0

At the moment, I am on Julia 1.1.1.

Do you have any idea why?

reworkhow commented 4 years ago

@fmorgante Can you update your JWAS now to the latest version to try?

fmorgante commented 4 years ago

@reworkhow It works now, thanks!

markcharder commented 4 years ago

Just a note to say thank you to the developers and I am greatly looking forward to trying out the software again :)