[R] uniroot
Petr Savicky
savicky at praha1.ff.cuni.cz
Sat Feb 5 08:38:42 CET 2011
On Fri, Feb 04, 2011 at 04:35:00AM -0800, dpender wrote:
>
> Hi,
>
> I am using the uniroot function in order to carry out a bivariate Monte
> Carlo simulation using the logistics model.
>
> I have defined the function as:
>
> BV.FV <- function(x,y,a,A)
> (((x^(-a^-1)+y^(-a^-1))^(a-1))*(y^(a-1/a))*(exp(-((1^(-a^-1)+y^(-a^-1))^a)+y^-1)))-A
>
> and the procedure is as follows:
>
> Randomly generate values of A~(0,1), y0 = -(lnA)^-1
>
> Where: A=Pr{X<xi|Y=yi-1} and a is the dependency between X and y (0.703)
>
> Use y0 to determine x where x = x(i) and y = y0(i-1) in the above equation.
>
> Unfortunately when the randomly defined A gets to approximately 0.46 the
> intervals I provide no longer have the opposite signs (both -ve).
>
> I have tried various different upper limits but I still can't seem to find a
> root.
>
> Does anyone have any suggestions?
Hi.
The expression for BV.FV() contains only one occurrence of x and it may be
separated from the equation. The solution for x may be expressed as
get.x <- function(y,a,A)
((A/(y^(a-1/a))/(exp(-((1^(-a^-1)+y^(-a^-1))^a)+y^-1)))^(1/(a-1))-y^(-a^-1))^(-a)
a <- 0.703
A <- 0.45
y <- - 1/log(A)
uniroot(BV.FV, c(1, 20), y=y, a=a, A=A)$root
[1] 3.477799
get.x(y, a, A)
[1] 3.477794
As pointed out in another reply, the equation does not always have
a solution. It seems to tend to infinity in the following
for (A in seq(0.45, 0.4677, length=20)) {
y <- - 1/log(A)
cat(get.x(y, a, A), uniroot(BV.FV, c(1, 1000), y=y, a=a, A=A)$root, "\n")
}
3.477794 3.477791
3.637566 3.637565
3.812851 3.812851
4.006213 4.006214
4.220838 4.220834
4.460738 4.460739
4.731048 4.731052
5.038453 5.038473
5.391845 5.391843
5.803329 5.803329
6.289873 6.28988
6.876084 6.876084
7.5992 7.599201
8.518665 8.518683
9.736212 9.736212
11.44329 11.44328
14.05345 14.05345
18.68155 18.68155
29.97659 29.9766
219.3156 219.3155
Hope this helps.
Petr Savicky.
More information about the R-help
mailing list