2

OBJECTIVE : 既知の行列 A と行列 C があります。C は未知の行列 BC の関数であり、A はほぼ等しいはずです。そのため、見つけたい未知の B を使用して AC を最小化しようとします。これには fminsearch を選択します。

1) 私のコスト関数では、ベクトル B を受け取り、それを使用してコスト関数を計算します。計算されたコスト関数は、ベクトルに変換してベクトルを返す行列です。

2) fminsearch を呼び出す/使用する直前に、B をベクトルとして定義し、B を初期化します (行列を定義し、ベクトルに変換します)。fminsearch を呼び出し、アルゴリズムの開始点として B を送信します。問題は滑らかではなく、不連続です。

5回繰り返した後、エラーが発生しました：

Iteration   Func-count     min f(x)         Procedure

 received B 0.1 0.1 0.1 0.1
 returning  0.163476621025508 0.293733517856535 0.127416423140963 -0.00271296077015587

 received B 0.105 0.1 0.1 0.1
 returning  0.16316485573574 0.293716885638579 0.127398072233325 -0.0027135923399549

 received B 0.1 0.105 0.1 0.1
 returning  0.163460015346026 0.293389431663331 0.127415792328966 -0.00273253871794082

 received B 0.1 0.1 0.105 0.1
 returning  0.163460015346026 0.293732946124453 0.127036780812566 -0.00273253871794082

 received B 0.1 0.1 0.1 0.105
 returning  0.1634760485781 0.293715744158374 0.127396812798301 -0.00313260869222808
     0            1        0.1634766                             
     1            5        0.1631649         initial simplex     

 received B 0.1025 0.1025 0.1025 0.1025
 returning  0.163303791321532 0.293543554409175 0.127206829827164 -0.00294271094775549
Error in this$fv[ive, 1] <- fv : 
  number of items to replace is not a multiple of replacement length
In addition: Warning messages:

コード フラグメントは次のとおりです。

require("neldermead")
require("fBasics")
  objective_f<-function(B){
  cat("\n received B",B)
  P_k_l_cts=matrix(0,nrow=regimes_optimum,ncol=regimes_optimum)
  B2=matrix(0,nrow=regimes_optimum,ncol=regimes_optimum)
  for(k in 1:(regimes_optimum*regimes_optimum))
  {
    m=ceiling(k/regimes_optimum)
    p=k%%regimes_optimum
    if(p==0)
    {p=regimes_optimum}
    B2[m,p]=B[k]
  }
  for(k in 1:regimes_optimum)
  {
    for(l in 1:regimes_optimum)
    {
      P_k_l_cts[k,l]=((omega%*%expm((B2-omega))*(1/lamda_cts[k]))[k,l])*(1/lamda_cts[k])
    }
  }
  P_k_l_cts_vec=vec(t(P_k_l_cts))
  P_k_l_d=matrix(0,nrow=regimes_optimum,ncol=regimes_optimum)
  for(k in 1:regimes_optimum)
  {
  for(l in 1:regimes_optimum)
  {
   if(k==l)
  {
    P_k_l_d[k,l]=(A_j_k_optimum[k,l])*(lamda_optimum[k])*(exp(-lamda_optimum[k]*    (1/lamda_optimum[k])))*(1/lamda_optimum[k])
   }
   else
   {
        P_k_l_d[k,l]=(lamda_optimum[k]/(lamda_optimum[k]-lamda_optimum[l]))*(A_j_k_optimum[k,l])*((exp(-lamda_optimum[l]*(1/lamda_optimum[k])))-(exp(-lamda_optimum[k]*(1/lamda_optimum[k]))))
      }
    }
  }

  P_k_l_d_vec=vec(t(P_k_l_d))
  cat("\n returning ",P_k_l_d_vec-P_k_l_cts_vec)
  cat("\n")
  return(P_k_l_d_vec-P_k_l_cts_vec)
}



    regime_cts_param<-function(A_j_k_cts,lamda_cts,p_cts,pi_cts,sigma_2_log_cts,regimes_optimum){      sigma_cts<<-matrix(0,nrow=regimes_optimum,ncol=1)
        for(j in 1:regimes_optimum)
     {
        sigma_cts[j]=sqrt(((1-p_cts[j])/lamda_cts[j])*exp(sigma_2_log_cts[j]))      
     }

  opt <- optimset(Display="iter",OutputFcn=outfun)
  output_fsolve=fminsearch(f=objective_f,x0=rep(0.1,regimes_optimum*regimes_optimum),opt)  
  return(list(sigma_cts=sigma_cts,B=output_fsolve$x,val=output_fsolve$fval))
}

内部の変数はグローバルであり、実際の値を提供しています:

regimes_optimum =2

sigma_2_log_optimum<-structure(c(-12.0947707925948, -12.0017466436778), .Dim = c(2L, 
1L))

lamda_optimum<-structure(c(2.1264772727783, 1.92731793847921), .Dim = c(2L, 
1L))

pi_optimum<-structure(c(0.792902540088065, 0.207097459911935), .Dim = c(2L, 
1L))

p_optimum<-structure(c(0.964696592041726, 0.393390784503164), .Dim = c(2L, 
1L))

A_j_k_optimum<-structure(c(0.613756531711541, 0.386243468288459, 0.779456839468863, 
0.220543160531137), .Dim = c(2L, 2L))

omega<-structure(c(2.1264772727783, 0, 0, 1.92731793847921), .Dim = c(2L, 
2L))

lamda_cts<-lamda_optimum
4

0 に答える 0