作者Dboy (废物)
看板R_Language
标题[问题] optim()
时间Fri Dec 6 10:47:26 2013
[软体熟悉度]:
入门(写过其他程式，只是对语法不熟悉)
[问题叙述]:
optim()侦测到的变数维度与预料中不同
各位板上先进好~
我现在是经济系的研究助理~
要帮老板用MitISEM这个package去做Baysian Estimation~
其中我的目标函数的程式码如下:
SChang <- function(theta, y, Z, w, t, dim, log=TRUE){
if (is.vector(theta))
theta <- matrix(theta, nrow = 1)
n <- length(y)
bx <- theta[,1]
bz <- matrix(theta[,2:(1+dim[1])],nrow=dim[1],ncol=1)
mu_x <- theta[,2+dim[1]]*matrix(rep(1,n),nrow=n, ncol=1)
del <- matrix(theta[,(3+dim[1]):(2+dim[1]+dim[2])],nrow=dim[2],ncol=1)
sig_v <- abs(theta[,3+dim[1]+dim[2]]) ## sig_v^2
sig_u <- abs(theta[,4+dim[1]+dim[2]]) ## sig_u^2
sig_epo <- abs(theta[,5+dim[1]+dim[2]]) ## sig_epo^2
sig_x <- abs(theta[,6+dim[1]+dim[2]]) ## sig_x^2
u <- matrix(theta[,(7+dim[1]+dim[2]):(6+dim[1]+dim[2]+n)],nrow=n, ncol=1)
x <- matrix(theta[,(7+dim[1]+dim[2]+n):(6+dim[1]+dim[2]+n+n)], nrow=n, ncol=1)
if (prod(u >0) > 0)
p <- -(n*log(sig_v^0.5)+n*log(sig_u^0.5)+n*log(sig_x^0.5)+n*log(sig_epo^0.5)+
(y+u-Z %*% bz-x*bx) %*% t(y+u-Z %*% bz-x*bx)/(2*sig_v)+
(w-x) %*% t(w-x)/(2*sig_epo)+
(x-mu_x) %*% t(x-mu_x)/(2*sig_x)+
(u-t %*% del) %*% t(u-t %*% del)/(2*sig_u)+
sum(log(pnorm(u/sig_u^0.5))))
else p <- -Inf
if (!log)
p <- exp(p)
as.vector(p)
}
其中y, Z, w, t是我的data~
theta在我的case中应该是个2014维的向量~
也是我要估计的函数的变数~
接下来我就把它丢到MitISEM中做估计:
app.Schang <- MitISEM(KERNEL=SChang, mu0=mu0, y=y, Z=Z, w=w, t=t, dim=dim)
结果显示：
错误在optim(par = mu0, fn = KERNEL, method = method, control = control, :
optim 内的目的函数长度被评估为 998001，而不是 1
是我目标函数写错吗?
998001与我预期的2014相去甚远~
我到底是哪里做错了?
[关键字]: optim()
--
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 140.112.4.183
※ 编辑: Dboy 来自: 140.112.4.183 (12/06 10:48)
1F：→ clickhere:p 回传的是 logL value, 不是vector of length(y). 12/06 11:28
2F：→ Dboy:这是MitISEM这个package的要求~它要求要有取log这个动作~所以 12/06 11:31
3F：→ Dboy:我才这样写~还是说我误会1楼的意思了? 12/06 11:31
4F：→ Dboy:而且optim本来就要求目标函数传回scalar而不是向量不是吗? 12/06 11:33
5F：→ Dboy:我发现我transpose写错地方......感谢各位~ 12/06 12:27