if (FALSE) {
####2-dim model
set.seed(123)
b <- c("-theta1*x1+theta2*sin(x2)+50","-theta3*x2+theta4*cos(x1)+25")
a <- matrix(c("4+theta5*sin(x1)^2","1","1","2+theta6*sin(x2)^2"),2,2)
true = list(theta1 = 0.5, theta2 = 5,theta3 = 0.3,
theta4 = 5, theta5 = 1, theta6 = 1)
lower = list(theta1=0.1,theta2=0.1,theta3=0,
theta4=0.1,theta5=0.1,theta6=0.1)
upper = list(theta1=1,theta2=10,theta3=0.9,
theta4=10,theta5=10,theta6=10)
start = list(theta1=runif(1),
theta2=rnorm(1),
theta3=rbeta(1,1,1),
theta4=rnorm(1),
theta5=rgamma(1,1,1),
theta6=rexp(1))
yuimamodel <- setModel(drift=b,diffusion=a,state.variable=c("x1", "x2"),solve.variable=c("x1","x2"))
yuimasamp <- setSampling(Terminal=50,n=50*100)
yuima <- setYuima(model = yuimamodel, sampling = yuimasamp)
yuima <- simulate(yuima, xinit = c(100,80),
true.parameter = true,sampling = yuimasamp)
prior <-
list(
theta1=list(measure.type="code",df="dunif(z,0,1)"),
theta2=list(measure.type="code",df="dnorm(z,0,1)"),
theta3=list(measure.type="code",df="dbeta(z,1,1)"),
theta4=list(measure.type="code",df="dgamma(z,1,1)"),
theta5=list(measure.type="code",df="dnorm(z,0,1)"),
theta6=list(measure.type="code",df="dnorm(z,0,1)")
)
mle <- qmle(yuima, start = start, lower = lower, upper = upper, method = "L-BFGS-B",rcpp=TRUE)
print(mle@coef)
set.seed(123)
bayes1 <- lseBayes(yuima, start=start, prior=prior,
method="mcmc",
mcmc=1000,lower = lower, upper = upper,algorithm = "randomwalk")
bayes1@coef
set.seed(123)
bayes2 <- lseBayes(yuima, start=start, prior=prior,
method="mcmc",
mcmc=1000,lower = lower, upper = upper,algorithm = "MpCN")
bayes2@coef
}
Run the code above in your browser using DataLab