我有一个函数,看起来像:g(x)= f(x)-a ^ b / f(x)^ b
g(x)-已知函数,提供的数据向量。
f(x)-隐藏过程。
a,b-此函数的参数。
从上面我们得到的关系:
f(x)= 逆(g(x))
我的目标是优化参数a和b,使 f (x) 尽可能接近正态分布。
我写了的代码:
g_fun <- function(x) {x - a^b/x^b}
inverse = function (f, lower = 0, upper = 2000) {
function (y) uniroot((function (x) f(x) - y), lower = lower, upper = upper)[1]
}
f_func = inverse(function(x) g_fun(x))
enter code here
# let's made up an example
# g(x) values are known
g <- c(-0.016339, 0.029646, -0.0255258, 0.003352, -0.053258, -0.018971, 0.005172,
0.067114, 0.026415, 0.051062)
# Calculate f(x) by using the inverse of g(x), when a=a0 and b=b0
for (i in 1:10) {
f[i] <- f_fun(g[i])
}
我有两个问题:
如何将参数 a 和 b 传递给函数?
如何执行此优化任务,这意味着找到 a 和 b,使得 f(x)近似正态分布。
不确定你如何能够产生 Q-Q 图,因为你提供的例子不起作用。你没有指定 a 和 b 的值,你正在定义 f_func,但调用 f_fun。无论如何这里是我对你的问题的回答:
如何将参数 a 和 b 传递给函数?-只需将它们作为参数传递给函数。
如何执行这个优化任务,意味着找到 a 和 b,使得 f(x)近似正态分布?-与任何优化任务的完成方式相同。定义一个成本函数,然后最小化它。
这是修订后的代码:我添加了 a 和 b 作为参数,删除了反函数并将其合并到 f_func 中,现在可以接受向量输入,因此不需要 for 循环。
g_fun <- function(x,a,b) {x - a^b/x^b}
f_func = function(y,a,b,lower = 0, upper = 2000){
sapply(y,function(z) { uniroot(function(x) g_fun(x,a,b) - z, lower = lower, upper = upper)$root})
}
# g(x) values are known
g <- c(-0.016339, 0.029646, -0.0255258, 0.003352, -0.053258, -0.018971, 0.005172,
0.067114, 0.026415, 0.051062)
f <- f_func(g,1,1) # using a = 1 and b = 1
#[1] 0.9918427 1.0149329 0.9873386 1.0016774 0.9737270 0.9905320 1.0025893
#[8] 1.0341199 1.0132947 1.0258569
f_func(g,2,10)
[1] 1.876408 1.880554 1.875578 1.878138 1.873094 1.876170 1.878304 1.884049
[9] 1.880256 1.882544
现在对于优化部分,这取决于你所说的 f (x) 将近似正态分布。如果你愿意,你可以比较 qq 线上的均方误差。既然你说近似,多近就足够了?你可以用 shapiro.test 继续搜索,直到你发现 p 值低于 0.05 (注意可能没有解决方案)
shapiro.test(f_func(g,1,2))$p
[1] 0.9484821
cost <- function(x,y) shapiro.test(f_func(g,x,y))$p
现在我们有了一个成本函数,我们如何去最小化它。有很多不同的方法来做数值优化。看看 optim 函数http://stat.ethz.ch/R-manual/R-patched/library/stats/html/optim.html。
optim(c(1,1),cost)
这最后一行不工作,但没有适当的数据和上下文,这是我可以去。
本站系公益性非盈利分享网址,本文来自用户投稿,不代表码文网立场,如若转载,请注明出处
评论列表(51条)