如何组织最终的 table 来模拟具有不同样本量的 r 中的最大似然?
how to organize the final table to simulate maximum likelihood in r with different sample sizes?
我使用 R 程序进行了模拟,但我无法 organize/make R 中的最终 table。
rep=100
n=20
# define starting values
mu<-100
sig<-3
theta=c(mu,sig) # store starting values
#Tables
#********
LHE=array(0, c(2, rep));
rownames(LHE)= c("MLE_mu", "MLE_sigma")
bias= array(0, c(2, rep));
rownames(bias)= c("bias_mu", "bias_sigma")
#Simulation {FOR LOOP}
#***********************
set.seed(1)
for(i in 1:rep){
myx <- rnorm(100, 100, 3)
loglikenorm<-function(x, myx) # always use x to hold parameter values
{
mu<-x[1]
sig<-x[2]
n<-length(myx)
loglike<- -n*log(sig*sqrt(2*pi))- sum((1/(2*sig^2))*(myx-mu)^2) # note
# use of sum
loglike<- -loglike
}
result<-nlm(loglikenorm, theta , myx=myx, hessian=TRUE, print.level=1) #ML estimation using nlm
mle<-result$estimate #extract and store mles
LHE[,i]= c(mle[1], mle[2])
bias[,i]= c(mle[1]-theta[1], mle[2]-theta[2])
} # end for i
L <-round(apply(LHE, 1, mean), 3) # MLE of all the applied iterations
bs <-round(apply(bias,1, mean),3) # bias of all the applied iterations
row<- c(L, bs); row
这将 运行 MLE 100 次。我想计算不同样本大小 (n=20,50,100) 和不同参数集的 MLE c(c(mu= 100, sigma=3), c(mu=80 , sigma=4))
我想要两件事,第一件事是如何 运行 代码来计算不同样本大小和不同参数集的 MLE。第二个如何使用 R 程序组织输出(如附图)。
The final table from R
我们将不胜感激。
我建议您构建一个用于 MLE 估计的函数,然后将其应用于不同的参数设置。这里的函数使用你的代码:
#Function
mymle <- function(n,mu,sig,rep)
{
#Set reps
rep=rep
# define starting values
mu<-mu
sig<-sig
theta=c(mu,sig) # store starting values
#Tables
LHE=array(0, c(2, rep));
rownames(LHE)= c("MLE_mu", "MLE_sigma")
#Bias
bias= array(0, c(2, rep));
rownames(bias)= c("bias_mu", "bias_sigma")
#Simulation
set.seed(1)
#Loop
for(i in 1:rep){
myx <- rnorm(rep, mu, sig)
loglikenorm<-function(x, myx) # always use x to hold parameter values
{
mu<-x[1]
sig<-x[2]
n<-length(myx)
loglike<- -n*log(sig*sqrt(2*pi))- sum((1/(2*sig^2))*(myx-mu)^2) # note
# use of sum
loglike<- -loglike
}
result<-nlm(loglikenorm, theta , myx=myx, hessian=TRUE, print.level=1) #ML estimation using nlm
mle<-result$estimate #extract and store mles
LHE[,i]= c(mle[1], mle[2])
bias[,i]= c(mle[1]-theta[1], mle[2]-theta[2])
} # end for i
#Format results
L <-round(apply(LHE, 1, mean), 3) # MLE of all the applied iterations
bs <-round(apply(bias,1, mean),3) # bias of all the applied iterations
row<- c(L, bs)
#Format a label
lab <- paste0('n= ',n,';',' mu= ',mu,';',' sig= ',sig)
row2 <- c(lab,row)
row2 <- as.data.frame(t(row2))
return(row2)
}
现在我们申请不同的参数设置:
#Example 1
ex1 <- mymle(n = 20,mu = 100,sig = 3,rep = 100)
ex2 <- mymle(n = 50,mu = 100,sig = 3,rep = 100)
ex3 <- mymle(n = 100,mu = 100,sig = 3,rep = 100)
#Example 2
ex4 <- mymle(n = 20,mu = 80,sig = 4,rep = 100)
ex5 <- mymle(n = 50,mu = 80,sig = 4,rep = 100)
ex6 <- mymle(n = 100,mu = 80,sig = 4,rep = 100)
最后,我们将所有结果绑定到接近您想要的输出:
#Bind all
df <- rbind(ex1,ex2,ex3,ex4,ex5,ex6)
输出:
V1 MLE_mu MLE_sigma bias_mu bias_sigma
1 n= 20; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
2 n= 50; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
3 n= 100; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
4 n= 20; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
5 n= 50; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
6 n= 100; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
我针对不同的参数组合分别测试了您的循环,结果相同。正如根据大数原理提醒的那样,随着rep
的增加,结果会趋向于真实值,所以如果将rep
改为1000,数值也会发生变化。
我使用 R 程序进行了模拟,但我无法 organize/make R 中的最终 table。
rep=100
n=20
# define starting values
mu<-100
sig<-3
theta=c(mu,sig) # store starting values
#Tables
#********
LHE=array(0, c(2, rep));
rownames(LHE)= c("MLE_mu", "MLE_sigma")
bias= array(0, c(2, rep));
rownames(bias)= c("bias_mu", "bias_sigma")
#Simulation {FOR LOOP}
#***********************
set.seed(1)
for(i in 1:rep){
myx <- rnorm(100, 100, 3)
loglikenorm<-function(x, myx) # always use x to hold parameter values
{
mu<-x[1]
sig<-x[2]
n<-length(myx)
loglike<- -n*log(sig*sqrt(2*pi))- sum((1/(2*sig^2))*(myx-mu)^2) # note
# use of sum
loglike<- -loglike
}
result<-nlm(loglikenorm, theta , myx=myx, hessian=TRUE, print.level=1) #ML estimation using nlm
mle<-result$estimate #extract and store mles
LHE[,i]= c(mle[1], mle[2])
bias[,i]= c(mle[1]-theta[1], mle[2]-theta[2])
} # end for i
L <-round(apply(LHE, 1, mean), 3) # MLE of all the applied iterations
bs <-round(apply(bias,1, mean),3) # bias of all the applied iterations
row<- c(L, bs); row
这将 运行 MLE 100 次。我想计算不同样本大小 (n=20,50,100) 和不同参数集的 MLE c(c(mu= 100, sigma=3), c(mu=80 , sigma=4))
我想要两件事,第一件事是如何 运行 代码来计算不同样本大小和不同参数集的 MLE。第二个如何使用 R 程序组织输出(如附图)。
The final table from R
我们将不胜感激。
我建议您构建一个用于 MLE 估计的函数,然后将其应用于不同的参数设置。这里的函数使用你的代码:
#Function
mymle <- function(n,mu,sig,rep)
{
#Set reps
rep=rep
# define starting values
mu<-mu
sig<-sig
theta=c(mu,sig) # store starting values
#Tables
LHE=array(0, c(2, rep));
rownames(LHE)= c("MLE_mu", "MLE_sigma")
#Bias
bias= array(0, c(2, rep));
rownames(bias)= c("bias_mu", "bias_sigma")
#Simulation
set.seed(1)
#Loop
for(i in 1:rep){
myx <- rnorm(rep, mu, sig)
loglikenorm<-function(x, myx) # always use x to hold parameter values
{
mu<-x[1]
sig<-x[2]
n<-length(myx)
loglike<- -n*log(sig*sqrt(2*pi))- sum((1/(2*sig^2))*(myx-mu)^2) # note
# use of sum
loglike<- -loglike
}
result<-nlm(loglikenorm, theta , myx=myx, hessian=TRUE, print.level=1) #ML estimation using nlm
mle<-result$estimate #extract and store mles
LHE[,i]= c(mle[1], mle[2])
bias[,i]= c(mle[1]-theta[1], mle[2]-theta[2])
} # end for i
#Format results
L <-round(apply(LHE, 1, mean), 3) # MLE of all the applied iterations
bs <-round(apply(bias,1, mean),3) # bias of all the applied iterations
row<- c(L, bs)
#Format a label
lab <- paste0('n= ',n,';',' mu= ',mu,';',' sig= ',sig)
row2 <- c(lab,row)
row2 <- as.data.frame(t(row2))
return(row2)
}
现在我们申请不同的参数设置:
#Example 1
ex1 <- mymle(n = 20,mu = 100,sig = 3,rep = 100)
ex2 <- mymle(n = 50,mu = 100,sig = 3,rep = 100)
ex3 <- mymle(n = 100,mu = 100,sig = 3,rep = 100)
#Example 2
ex4 <- mymle(n = 20,mu = 80,sig = 4,rep = 100)
ex5 <- mymle(n = 50,mu = 80,sig = 4,rep = 100)
ex6 <- mymle(n = 100,mu = 80,sig = 4,rep = 100)
最后,我们将所有结果绑定到接近您想要的输出:
#Bind all
df <- rbind(ex1,ex2,ex3,ex4,ex5,ex6)
输出:
V1 MLE_mu MLE_sigma bias_mu bias_sigma
1 n= 20; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
2 n= 50; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
3 n= 100; mu= 100; sig= 3 99.98 3.015 -0.02 0.015
4 n= 20; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
5 n= 50; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
6 n= 100; mu= 80; sig= 4 79.974 4.02 -0.026 0.02
我针对不同的参数组合分别测试了您的循环,结果相同。正如根据大数原理提醒的那样,随着rep
的增加,结果会趋向于真实值,所以如果将rep
改为1000,数值也会发生变化。