/
random.jl
86 lines (80 loc) · 3.39 KB
/
random.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
"""
rand(rng::AbstractRNG, lmm::LMM{T}) where T
Generate random responce vector for fitted 'lmm' model.
"""
function rand(rng::AbstractRNG, lmm::LMM)
if !lmm.result.fit error("Model not fitted!") end
rand(rng, lmm, lmm.result.theta, lmm.result.beta)
end
function rand!(rng::AbstractRNG, v::AbstractVector, lmm::LMM)
if !lmm.result.fit error("Model not fitted!") end
rand!(rng, v, lmm, lmm.result.theta, lmm.result.beta)
end
"""
rand!(v::AbstractVector, lmm::LMM) = rand!(default_rng(), v, lmm, lmm.result.theta, lmm.result.beta)
Generate random responce vector for fitted 'lmm' model, store results in `v`.
"""
rand!(v::AbstractVector, lmm::LMM) = rand!(default_rng(), v, lmm, lmm.result.theta, lmm.result.beta)
"""
rand(rng::AbstractRNG, lmm::LMM{T}; theta) where T
Generate random responce vector 'lmm' model, theta covariance vector, and zero means.
"""
function rand(rng::AbstractRNG, lmm::LMM{T}, theta::AbstractVector) where T
v = Vector{T}(undef, nobs(lmm))
rand!(rng, v, lmm, theta)
end
"""
rand!(rng::AbstractRNG, lmm::LMM{T}; theta) where T
Generate random responce vector 'lmm' model, theta covariance vector, and zero means, store results in `v`.
"""
function rand!(rng::AbstractRNG, v::AbstractVector, lmm::LMM{T}, theta::AbstractVector) where T
n = length(lmm.covstr.vcovblock)
if length(v) != nobs(lmm) error("Wrong v length!") end
tv = Vector{T}(undef, lmm.maxvcbl)
gvec = gmatvec(theta, lmm.covstr)
rtheta = lmm.covstr.tr[lmm.covstr.rn + 1:end]
for i = 1:n
q = length(lmm.covstr.vcovblock[i])
V = zeros(q, q)
vmatrix!(V, gvec, theta, rtheta, lmm, i)
if length(tv) != q resize!(tv, q) end
Distributions.rand!(rng, MvNormal(Symmetric(V)), tv)
v[lmm.covstr.vcovblock[i]] .= tv
end
v
end
rand(lmm::LMM, theta::AbstractVector) = rand(default_rng(), lmm, theta)
rand!(v::AbstractVector, lmm::LMM, theta::AbstractVector) = rand!(default_rng(), v, lmm, theta)
"""
rand(rng::AbstractRNG, lmm::LMM{T}; theta, beta) where T
Generate random responce vector 'lmm' model, theta covariance vector and mean's vector.
"""
function rand(rng::AbstractRNG, lmm::LMM{T}, theta::AbstractVector, beta::AbstractVector) where T
v = Vector{T}(undef, nobs(lmm))
rand!(rng, v, lmm, theta, beta)
end
function rand!(rng::AbstractRNG, v::AbstractVector, lmm::LMM{T}, theta::AbstractVector, beta::AbstractVector) where T
if length(beta) != size(lmm.data.xv, 2) error("Wrong beta length!") end
if length(theta) != lmm.covstr.tl error("Wrong theta length!") end
if length(v) != nobs(lmm) error("Wrong v length!") end
n = length(lmm.covstr.vcovblock)
tv = Vector{T}(undef, lmm.maxvcbl)
m = Vector{T}(undef, lmm.maxvcbl)
gvec = gmatvec(theta, lmm.covstr)
rtheta = lmm.covstr.tr[lmm.covstr.rn + 1:end]
for i = 1:n
q = length(lmm.covstr.vcovblock[i])
if length(tv) != q
resize!(tv, q)
resize!(m, q)
end
mul!(m, lmm.dv.xv[i], beta)
V = zeros(q, q)
vmatrix!(V, gvec, theta, rtheta, lmm, i)
Distributions.rand!(rng, MvNormal(m, Symmetric(V)), tv)
v[lmm.covstr.vcovblock[i]] .= tv
end
v
end
rand(lmm::LMM, theta::AbstractVector, beta::AbstractVector) = rand(default_rng(), lmm, theta, beta)
rand!(v::AbstractVector, lmm::LMM, theta::AbstractVector, beta::AbstractVector) = rand!(default_rng(), v, lmm, theta, beta)