# Calculating probability of displacement using two CDFs

My knowledge of stats is fairly basic, so you please bear with me!

I'm trying to calculate the CDF for the vertical displacement of a (light, small) object floating in a wave tank. I have experimental data which suggests that the displacement can be approximated to a normal distribution, however I am also required to provide an analytical solution, as the experimental data is not reliable enough.

My approach so far has been to use the Rayleigh Distribution to calculate the probability of any given wave being greater than an amplitude (h). $$P(H>h) = e^{-2\left(\frac{h^2}{h_s^2}\right)}$$ where $$h_s = 2\sigma$$

Separately, I have approximated each individual wave to a sin wave of constant frequency. I have calculated the CDF of a sin wave (of amplitude H) to be: $$P(X>x) = \frac{1}{2}-\frac{sin^{-1}(\frac{x}{H})}{\pi}$$ (where x is the vertical displacement of the object)

The part I do not understand, is how I combine these two to get the probability of exceeding a certain displacement for any wave. Instinctively, I suspect I am supposed to multiply the two and integrate for h between x and inf, but so far I haven't found an answer which is anywhere near what I am expecting, so I may be way off!

• How are $H$ and $X$ related? Without this information, you cannot specify the joint distribution. Jun 2, 2016 at 14:10
• For an individual wave, $H$ is the amplitude of that wave (I'll update my question as I think I put $h$ by mistake). The displacement at a point on a wave ($X$) obviously depends on the amplitude of that wave ($H$), $X=Hsin(\theta)$, but it is also dependent on the point on the wave $(\theta)$. Is that what you were looking for? Jun 2, 2016 at 14:22
• This begs for another question then: what is the distribution of $\theta$ and in particular can it be assumed to be independent of $X$? Jun 2, 2016 at 14:53
• The waves are assumed to be passing the object with constant frequency WRT time, so I suppose $\theta$ should really be $\omega.t$, where $\omega$ is the period (constant), and $t$ is time (so independent of $X$) Jun 2, 2016 at 14:58

Try to define the probability that you want to obtain in mathematical terms. If you do so the following becomes just calculus.

The probability you want to obtain can be seen as an expectation of $I\left[x > x_0\right]$ over joint distribution of $x$ and $h$ where $I$ is an indicator function that is equal to $1$ when $x > x_0$ and to $0$ otherwise: $$P(x > x_0) = \mathop{\mathbb E} _{x, h} I\left[x > x_0\right] = \int dP(x, h)\ I\left[x > x_0\right] = \iint dx\ dh\ p(x, h)\ \theta(x - x_0).$$

In this expression $p(x,h)$ is a probability density of joint distribution of $x$ and $h$, and $\theta(x)$ is a Heaviside Step function.

Your probability $P(X > x)$ is actually a conditional probability for given value of $h$ and its derivative is a probability density of conditional probability of $x$ for given $h$ $p(x \mid h).$ So using definition of conditional probability you can get that $p(x, h)\ dx\ dh = p(x \mid h)\ dx \cdot p(h)\ dh.$

So now you can obtain the probability density functions $$p(x \mid h) = \frac{d}{dx} P(X \leq x \mid h) = \frac{d}{dx}\left(1 - P(X > x \mid h) \right) = -\frac{\partial}{\partial x}P(X > x \mid h),$$ $$p(h) = \frac{d}{dh} P(H \leq h) = -\frac{\partial}{\partial h} P(H > h)$$ and substitute them to the integral: $$P(x > x_0) = \iint dx\ dh\ p(x\mid h)\ p(h)\ \theta(x - x_0)= \int\limits_0^{\infty} dh\ p(h) \int\limits_{x_0}^{\infty} dx\ p(x\mid h) = \\-\int\limits_0^{\infty} dh\ P(X > x_0 \mid h) \frac{\partial}{\partial h} P(H > h).$$

Actually in this particular case there is a more simple way to use conditional probabilities directly: $$P(X > x) = \mathop{\mathbb E}_h P(X > x \mid h) = \int P(X > x\mid h) \ dP(h),$$ but I find the joint distribution/indicator function approach to be more powerful because it allows to infer probability distribution of more complicated relations (e.g. if probability of wave energy $x^2$ is larger than given function of $h$) or introduce new probabilistic variables to the model (e.g. instead of single value of $\sigma$ use random variable sigma from distribution obtained in result of bayesian estimation of $\sigma$ from experimental data) using a single generic framework.

• That looks like exactly what I need, it'll take me a little while to work through but once done I'll mark as answered. Thanks a lot!! Jun 3, 2016 at 7:42
• So from my understanding of this, I get something like this: $$\\ P(X>x|h)=\frac{1}{2}-\frac{1}{\pi} sin^{-1}\left(\frac{x}{h} \right ) \\ \\ P(H>h)=e^{-2\frac{h^2}{h_s^2}} \\ \\ \frac{\partial}{\partial h}P(H>h)=\frac{4h}{h_s^2}e^{-2\frac{h^2}{h_s^2}} \\ \\ P(x>x_0)=-\int_0^\infty\left [ \frac{1}{2}-\frac{1}{\pi} sin^{-1}\left(\frac{x_0}{h} \right ) \right ]\left[ \frac{4h}{h_s^2}e^{-2\frac{h^2}{h_s^2}}\right ] dh$$. Does this look about right? Then I just need to work out the integral! Jun 3, 2016 at 8:57
• Note that you should integrate from $x_0$ to $\infty$, otherwise $\sin^{-1}$ is not defined if $h < x_0.$ Btw it's unlikely that this integral is expressible in elementary functions, but you can investigate numerically how does it behave making a substitution $h = h_s \cdot \xi$ where $\xi$ is a dimensionless variable and considering it as a function of a single parameter $x_0 / h_s.$ Jun 3, 2016 at 9:59
• Yeah you're right I can't edit the comment but I will integrate from $x_0$ to $\infty$. I agree that it isn't expressible in elementary functions. I think that actually makes sense, as I am aiming to get to something close to a normal distribution CDF, which also isn't. I'll have a go at a numerical solution, and post anything I get here. Either way, I've accepted your answer as it was exactly what I needed. Thanks so much! Jun 3, 2016 at 10:01
• If you replace $h$ to $\xi h_s$ than you can obtain some asymptotics using Tayler series expansion by $x_0/h_s$ parameter if you have reasons to think that this parameter is small (or by $h_s/x_0$ if you think that it is large) which could be close to normal. Also I recommend you to work with probability $P(x \leq x_0) = 1 - P(x > x_0)$ because it seems that calculations are simpler with it. Jun 3, 2016 at 10:10