作者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