有关R

问题描述

我对R中的贝叶斯最优设计有疑问。

我不确定我是否正确地遵循了内容

所以,对于贝叶斯最优性的基本思想,我认为

选择可最大化以下功能的设计计划:

enter image description here

然后我通过“贝叶斯最优设计用于具有部分轮廓的离散选择实验”修改了论文的功能

enter image description here

问题是我一直对每个集成部分感到困惑。

数据如下所示,我进行了普通的D最优设计(不是贝叶斯设计)

带有以下代码

顺便说一下,数据显示如下。

enter image description here

#Parameter estimation for the models

total.w<-function(param,data){
a1<-param[1]
b1<-param[2]
c1<-param[3]
value<-c()
for (i in 1:dim(data)[1]){
value[i]<-data[i,]$failure*log(1-exp(-(data[i,]$time/exp(a1+b1*data[i,]$stress))^c1))+(data[i,]$total-data[i,]$failure)*((-(data[i,]$stress))^c1))}
return(-sum(value))}
opt.ans4<-optim(par=ini.w,total.w,data=dataF2,hessian=TRUE)


#Making Fisher information matrix for optimal design criterion


  
a1<-opt.ans4$par[1]

b1<-opt.ans4$par[2]



c1<-opt.ans4$par[3]

id<-matrix(c("a1","a1","b1","c1","c1"),ncol=2,byrow=TRUE)

total.w<-function(tsample,fail,con,labT,a1,b1,c1){
fail*log(1-exp(-(labT/exp(a1+b1*con))^c1))+(tsample-fail)*((-(labT/exp(a1+b1*con))^c1))
}


ele1.w<- Deriv(Deriv(total.w,id[1,1]),2])
ele2.w<- Deriv(Deriv(total.w,id[2,2])
ele3.w<- Deriv(Deriv(total.w,id[3,2])
ele4.w<- Deriv(Deriv(total.w,id[4,2])
ele5.w<- Deriv(Deriv(total.w,id[5,2])
ele6.w<- Deriv(Deriv(total.w,id[6,2])
ele7.w<- Deriv(Deriv(total.w,id[7,2])
ele8.w<- Deriv(Deriv(total.w,id[8,2])
ele9.w<- Deriv(Deriv(total.w,id[9,2])

w1<- sum(-ele1.w(dataF2$total,dataF2$failure,dataF2$stress,dataF2$time,c1))
w2<- sum(-ele2.w(dataF2$total,c1))
w3<- sum(-ele3.w(dataF2$total,c1))
w4<- sum(-ele4.w(dataF2$total,c1))
w5<- sum(-ele5.w(dataF2$total,c1))
w6<- sum(-ele6.w(dataF2$total,c1))
w7<- sum(-ele7.w(dataF2$total,c1))
w8<- sum(-ele8.w(dataF2$total,c1))
w9<- sum(-ele9.w(dataF2$total,c1))
#
amat3<-matrix(c(w1,w2,w3,w4,w5,w6,w7,w8,w9),ncol=3)
opt.ans4$hessian


##Function to generate failure rate for experimental design

fail.func<-function(tsample,c1){
tsample*(1-exp(-(labT/exp(a1+b1*con))^c1))
}

#Example 
param<-c(8,8,10,10)
round(fail.func(dataF2$total,c1)) 

###


Test.w<-function(param){

samp<-round(c(param[1],param[2],param[3]))
freq<-round(5*c(param[4],param[5],param[6]))

#a1<- 5.2265
#b1<- -0.0351
#c1<- 1.9547
a1<-4.150
b1<- -1.384
c1<-1.2142

temp<-c(0.4286,0.7143,1.000)
w1.e<-c()
w2.e<-c()
w3.e<-c()
w4.e<-c()
w5.e<-c()
w6.e<-c()
w7.e<-c()
w8.e<-c()
w9.e<-c()
failure.w<-c()
for (i in 1:length(temp)){
failure.w[i]<-fail.func(param[i],temp[i],freq[i],c1)
}

for (i in 1:length(temp))
{
w1.e[i]<- (-ele1.w(samp[i],(failure.w[i]),c1))
w2.e[i]<- (-ele2.w(samp[i],c1))
w3.e[i]<- (-ele3.w(samp[i],c1))
w4.e[i]<- (-ele4.w(samp[i],c1))
w5.e[i]<- (-ele5.w(samp[i],c1))
w6.e[i]<- (-ele6.w(samp[i],c1))
w7.e[i]<- (-ele7.w(samp[i],c1))
w8.e[i]<- (-ele8.w(samp[i],c1))
w9.e[i]<- (-ele9.w(samp[i],c1))
}

fish<-matrix(c(sum(w1.e),sum(w2.e),sum(w3.e),sum(w4.e),sum(w5.e),sum(w6.e),sum(w7.e),sum(w8.e),sum(w9.e)),ncol=3)
value<- det(fish)
#print(fish)
#print(value)
return(value)
}

我的问题是

  1. 对于决定渔民信息的人,如何将功能修改为?我是否只是像下面的代码那样将先验分布函数放入对数似然函数中,并从中提取出更多信息?

     total.w<-function(tsample,c1){
    

    失败 log(1-exp(-(labT / exp(a1 + b1 con))^ c1))+(tsample-fail)((-(labT / exp( a1 + b1 con))^ c1)) +优先功能1+优先功能2+优先功能3}

  2. 第二个条件是为D最优生成设计空间。我对每个样本sie使用迭代,并为每次检查使用时间。是否与y的贝叶斯D最优的外部积分相同?还是我需要为此应用mcmc函数

解决方法

暂无找到可以解决该程序问题的有效方法,小编努力寻找整理中!

如果你已经找到好的解决方法,欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@)