如何在R中最大化函數內部使用矩陣乘法進行優化工作

我試圖最大化維數爲2x2的矩陣參數的可能性。似然函數需要傳遞幾個固定的矩陣參數，這個可能性也是一個函數。數據表示爲Y，協方差矩陣Sigma.star（我作爲下三角矩陣傳遞）是計算所必需的，但我希望保留這些數據並對其進行優化，在我的代碼試圖優化A如何在R中最大化函數內部使用矩陣乘法進行優化工作

我的問題是，似乎優化似乎是錯誤的事實，它優化了我用於矩陣代數的對象內部的東西。有沒有辦法讓它工作，而不需要編程每一個小小的計算呢？

特定的錯誤是：

Error in diag(1, nrow = (m^2)) - A %x% A : non-conformable arrays

但克羅內克一個是m^2 XM^2矩陣就像身份......

代碼：

library(MCMCpack) 
library(mvtnorm) 
set.seed(1000) 


Likelihood.orig<-function(A, Y, Sigma.star){ 
Sigma<-xpnd(Sigma.star) 
n<-nrow(Y) 
if(is.vector(A)==TRUE){ 
    A<-as.matrix(A, nrow=nrow(Sigma), ncol=ncol(Sigma)) 
} 
m<-nrow(A) 
V<-matrix(solve(diag(1, nrow=(m^2))-A%x%A)%*%as.vector(Sigma), nrow=m, ncol=m) 
temp1<- (-.5)*log(abs(det(V))) 
temp2<- (-(n-1)/2)*log(abs(det(Sigma))) 
temp3<- t(Y[,1, drop=FALSE]) %*% (solve(V)) %*% Y[,1, drop=FALSE] 
terms<- numeric(n-1) 
for(i in 2:n){ 
    terms[i-1]<- t(Y[,i, drop=FALSE] - A %*%Y[,i-1, drop=FALSE]) %*% (solve(Sigma)) %*% (Y[,i] - A %*%Y[,i-1]) 
} 
return(temp1+temp2-.5*(temp3+sum(terms))) 
} 




Generate.Y<-function(n, A, Sigma){ 
m<-nrow(A) 
Y<-matrix(0, nrow=m, ncol=n) 
V<-matrix(solve(diag(1, nrow=m^2)-A%x%A)%*%as.vector(Sigma), nrow=m, ncol=m) 
Y[,1]<-rmvnorm(1, numeric(nrow(A)), V) 
for(i in 2:n){ 
    Y[,i]<-A%*%Y[,i-1, drop=FALSE]+t(rmvnorm(1, mean = numeric(m), sigma = Sigma)) 
    } 
return(Y) 
} 


n<-500 
A.true<-matrix(c(.8, .3, 0, .5), nrow=2, ncol=2) 
Sigma<-matrix(c(1, 0, 0, .5), nrow=2, ncol=2) 
Y<-matrix(0, nrow=2, ncol=n) 
Y<-Generate.Y(n, A.true, Sigma) 
m=nrow(Y) 
lower.Sigma<-vech(Sigma) 



optim(par=c(1, 0, 0, 1), fn=Likelihood.orig, method="Nelder-Mead", 
    control=list(maxit=500, fnscale=-1), Sigma.star=lower.Sigma, Y=Y)

來源

2012-11-14 Peter

你的方法是正確的，即，使01優化一個向量，並且只將該向量轉換爲您嘗試最大化的函數內的矩陣。

但是，您需要使用matrix而不是as.matrix來創建該矩陣。見的區別：

as.matrix(1:4, nrow=2, ncol=2) # wrong tool 
#  [,1] 
# [1,] 1 
# [2,] 2 
# [3,] 3 
# [4,] 4

和

matrix(1:4, nrow=2, ncol=2) 
#  [,1] [,2] 
# [1,] 1 3 
# [2,] 2 4

對於這種類型的問題，我會強烈建議你學習將R調試工具（browser，debug，debugonce等）。例子見General suggestions for debugging in R。

來源

2012-11-15 00:02:14 flodel

如何在R中最大化函數內部使用矩陣乘法進行優化工作

回答

相關問題