Quantcast

uniroot

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

uniroot

dpender
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?

Thanks,

Doug
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: uniroot

Berend Hasselman
<quote author="dpender">
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).
</quote>

Try

curve(BV.FV(x,y=y0,a=.703,A=.7),from=1,to=10)

for different values of a, A and y0.
I have a feeling that your function or something else is not quite correct.

Berend
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: uniroot

Petr Savicky-2
In reply to this post by dpender
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.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: uniroot

dpender
Thanks for your advice.  There was an error in the equation that is was copying.

Doug
Loading...