Thanks for your reply!
I read the documentation but still not sure if this would work in my case (sorry I’m still new to Julia!). Let me describe what I am trying to do.
I defined the following function
function myfunc(δ, θ2, y, c, d, ϵ1, ϵ2)
N = size(y, 1)
βY = θ2[1:5]
γ0 = θ2[6]
γ1 = θ2[7:11]
τ = θ2[12]
κ = θ2[13]
s = d + 1/sqrt(τ) * ϵ1
p = p0./(p0 .+ (ones(N) - p0) .* exp.((0.5*ones(N) - s).*τ))
Eu = δ .+ y*βY - (γ0 .+ y*γ1).*(ones(N) - c).*p .+ ϵ2 .- (log(κ) - 1)
b = isequal.(sign.(Eu), 1.0)
return b
end
The input are: δ (scaler), θ2 (13-by-1 vector), y (N-by-5 matrix), c, d, ϵ1, ϵ2 (N-by-1 vector). This function returns a N-by-1 vector, and N is around 1000. What I really want is a vector whose elements are the expectation of b’s elements over ϵ1, ϵ2, which are standard normal variables and mutually independent.
Is it possible to do the integration within the function, so instead of having ϵ1, ϵ2 as inputs, having the function directly return the calculated expectations?