作者f496328mm (为什麽会流泪)
看板R_Language
标题[问题] MLE 做 联立方程式 估计 (已解决)
时间Wed Sep 7 15:56:29 2016
联立方程式 or system equation
用MLE去做 我目前只会做一条equation的MLE
我目前上网找的只能做到一条equation的MLE估计
目前code如下:
#==========================================
library(stats4)
set.seed(1001)
N <- 100
x <- runif(N)
y <- 5 * x + 3 + rnorm(N)
LL <- function(beta0, beta1, mu, sigma) {
# Find residuals
#
R = y - x * beta1 - beta0
#
# Calculate the likelihood for the residuals (with mu and sigma as parameters)
#
R = suppressWarnings(dnorm(R, mu, sigma))
#
# Sum the log likelihoods for all of the data points
#
-sum(log(R))
}
fit <- mle(LL, start = list(beta0 = 3, beta1 = 1, mu = 0, sigma=1))
summary(fit)
#==========================================
R = y - x * beta1 - beta0
这是我的一条equation
我目前只会做到一条
如果要做一条以上，要如何做?
#=========================================
另外我也有找到这个package "systemfit"
它可以一次做两个以上的回归
不过 我的理解是 他是分开做的
因为我分别用lm去做 估计出来的参数也一样
所以没什麽意义 可能是我理解错误
不知道这个函数的用意是什麽
谢谢
--
※ 发信站: 批踢踢实业坊(ptt.cc), 来自: 114.42.98.232
※ 文章网址: https://webptt.com/cn.aspx?n=bbs/R_Language/M.1473234992.A.2CE.html
1F：推 ksherry: Google Full information maximum likelihood 09/07 17:23
2F：→ ksherry: 联立方程组 你还是要自己先求出均衡解 09/07 17:23
3F：推 qhair: 我猜你要的是BBsolve 09/07 17:39
不是BBsolve喔 刚刚看了一下 不过还是谢谢
因为是估计值 所以至少要有Std. Error p-value
我正在试试看分开用nls估计 看差异大不大
只要还是想一起估计 因为互相为限制式 可是我不会多条equation估计><
※ 编辑: f496328mm (114.42.98.232), 09/07/2016 19:05:26
自己的问题自己回
研究出来啦 烦了我好几个礼拜
http://imgur.com/xUatDl5
先做模拟预设要估计的beta 两条equation
最後估计出近似原始设的beta
code如下:
library(stats4)
set.seed(1001)
N <- 100
x1 <- rnorm(N,10,1)
x2 <- rnorm(N,5,1)
beta0=3
beta1=5
beta2=50
y1 <- beta1 * x1 + beta0 + rnorm(N)
y2 <- beta2 * x2 + beta1*x1 + beta0 + rnorm(N)
LL <- function(beta0, beta1,beta2, mu, sigma) {
R1 = y1 - x1 * beta1 - beta0
R2 = y2 - x2 * beta2 - x1*beta1 - beta0
#
R1 = suppressWarnings(dnorm(R1, mu, sigma, log = TRUE))
R2 = suppressWarnings(dnorm(R2, mu, sigma, log = TRUE))
#
-sum(R1,R2)
}
fit <- mle2(LL, start = list( beta0 = 3, beta1 = 1,
beta2 = 1, mu = 0, sigma = 1))
summary(fit)
※ 编辑: f496328mm (114.42.98.232), 09/07/2016 19:42:41
※ 编辑: f496328mm (114.42.98.232), 09/07/2016 19:42:54
4F：→ celestialgod: 发文格式... 09/07 20:19