5 votos

Región de parámetro para la existencia de soluciones de la ecuación

¿Cómo funciona una parcela en el set de parámetros región de ecuaciones usando el r?. Mi problema es que tengo una ecuación con 2 variables definidas por $g(x,a,b) = 0$ y quiero graficar el conjunto de $(a,b)$ cuando la solución existe.

El ejemplo concreto que me interesa es la siguiente:

$$(F(x) - a)f^2(x) + \left(\frac{F^2(x)}{2} - a F(x) + \left(\frac{a^2}{2} + b\right)\right)f^\prime(x)=0.$$

where $F$ is the cumulative function of a normal distribution of mean=0 and sd=1, $f$ the density function of the normal, and $f^\prime$ is the first derivative of $f$.

My goal is to construct the region (surface) $S$ where

$$S = \{(a,b)\in [0,1]^2\ |\text{ a solution }x\text{ of the equation exists and belongs to }[0,1]\}. $$

E<-function(a,u) {  
((a-u)*(dnorm(q(u),m,s, log = FALSE))/(qnorm(u)))-((qnorm(u))^2)/2 +a*qnorm(u)-(a^2)/2 
}
a <- u <- seq(0,1,0.01)
z <- outer(a,u,E)
b<-seq(0,1,0.001)
contour( x=a, y=a, z=z,levels=b, las=1, drawlabels=FALSE, lwd=3,xlab="a", ylab="x")

Please help me here. I hope this can teach me more. Thanks in advance.


Thanks, this is exactly what I was looking for. But I have 2 remarks or comments : Remark 1 : I made a mistake $u=F(x)\in[0,1]$ but $x = N(0,1)$ so $x$ is not restricted in $[0,1]$. Por lo tanto, es posible proponer que el código que se utiliza en Mathematica tal que yo pueda tratar de adaptarla.

Observación 2 : me hizo también una simulación usando Mathematica y me sale la imagen en 3D de mi ecuación, pero no fue posible obtener la proyección para obtener S. Este es el código que he usado :

f[u_] := Quantile[NormalDistribution[0, 1], u]
g[u_]:=PDF[NormalDistribution[0,1],f[u]]
h[u_, a_, e_] := (u - a)*((g[u])^2) + ((u^2)/2 - a*u + (a^2)/2 + e)*f[u]*g[u]
ContourPlot3D[h[u, a, e] == 0, {a, 0, 1}, {u, 0, 1}, {e, 0, 1}]

Luego me sale esta imagen :

enter image description here

Entonces, ¿hay forma de obtener la proyección de esta imagen.

5voto

jldugger Puntos 7490

Para abordar la pregunta general que considere la posibilidad de usar una herramienta que se adapta bien a este tipo de cálculos y visualizaciones, tales como Mathematica. (Esto fue utilizado para la trama de las dos primeras y dos últimas cifras, a continuación.)

Esta pregunta en particular es susceptible de ulterior análisis, que permite a R a mostrar $S$:$x\in [0,1]$, podemos determinar el subconjunto de $S$ determina (que es una curva). Por la elección de un visualmente densa colección de tales $x$, la colección de estas curvas limns toda la región de $S$.


Comenzar el análisis por la reescritura de la definición de la ecuación en la forma

$$b = -\frac{1}{2} a^2 + u(x) a + v(x)$$

where (assuming $x\ne 0$)

$$u(x) = f(x)^2/f^\prime(x) + F(x)$$ and $$v(x) = -\left(F(x) f(x)^2/f^\prime(x) + \frac{1}{2}F(x)^2\right).$$

In the $(a,b)$ plane these equations describe similar parabolae having their vertexes at $(u(x), v(x) + u(x)^2/2).$

Figure 1

In this figure some of the parabolae are drawn in gray. The locus of their vertexes is traced by the thick red curve. The region $[0,1]\los tiempos de [0,1]$ is shown as a gray square. $S$ is the portion of the gray square overlapped by the shaded blue region.

The region $S$ comprises most of the left half of the unit square ($\le 1/2$), less a small region at the bottom, together with pieces of some parabolae at the bottom. The next figure examines that lower region in more detail.

Figure 2

The piece missing from $S$ in the left half of the unit square lies below the parabola $b = -\frac{1}{2}a^2 + u(1) + v(1)\aprox -\frac{1}{2}a^2 + 0.599374 un - 0.15035.$ The part of $S$ in the right half of the unit square (where the curves seem to overlap) is bounded by the envelope of these parabolae; it is difficult to derive any simple formula for its boundary.


R code

Begin with a function f to compute values along a parabola given by $x \ne 0$:

f <- function(a, x) {
  F <- pnorm(x)
  f0 <- dnorm(x)
  f1 <- -x * exp(-x^2 / 2) / sqrt(2 * pi)
  return(-1/2 * a^2 + (f0^2/f1 + F)*a - (F*f0^2/f1 + F^2/2))
}

Because interest lies in $0\le\le 1$, a typically will be a sequence of numbers in this range; f returns the ordinates of the parabola lying above this sequence.

Use this iteratively to draw the parabolae. The vertical line $a=\frac{1}{2}$ (corresponding to $x=0$) needs to be drawn separately because it is not the graph of a function of $$.

n.mesh <- 64
mesh <- seq(0, 1, length.out=n.mesh)
colors <- hsv(mesh, .8, .9, 2/3)
plot(c(0,1), c(0,1), type="n", xlab="a", ylab="b")
rect(0,0,1,1, col=gray(0.96))

for(i in n.mesh:2) {
  x0 <- mesh[i]
  curve(f(x, x0), add=TRUE, col=colors[i])
}
abline(v=1/2, xlim=c(0,1), ylim=c(0,1), col=colors[1])

Figure 3


Appendix: Brute Force Code

In Mathematica a brute-force (but somewhat efficient) way to find a subset of $S$, when the function $x\g(x,a,b)$ is continuous for all $(a,b)\en[0,1]^2$, checks whether this function changes sign between $x=0$ and $x=1$. Begin by defining $g$:

fF[x_] := CDF[NormalDistribution[], x];
f[x_] := PDF[NormalDistribution[]][x];
f1[x_] := Evaluate[D[f[y], y] /. y -> x];
g[x_, a_, b_] := (fF[x] - a) f[x]^2 + (fF[x]^2/2 - a fF[x] + (a^2/2 + b)) f1[x];

Here is the solution (which plots quickly):

RegionPlot[g[0, a, b] g[1, a, b] < 0 , {a, 0, 1}, {b, 0, 1}]

Figure 3

In general some analysis is required to assure that this covers all of $S$. That analysis can be assisted by plotting contours of $(a,b)\g(x,a,b)$ for various values of $x$ de interés:

ContourPlot[Evaluate@Table[g[x, a, b] == 0, {x, Range[0, 1, 1/64]}], {a, 0, 1}, {b, 0, 1}]

Figure 4

En la medida en la que puede confiar R para producir precisa el contorno de las parcelas, una solución similar para el mismo.

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X