# Fast integration of 2D distribution across lines parallel to y-axis

I'm struggling with a small data set and a slow calculation. I have hundreds of small 2D data arrays and need to integrate across several lines parallel to the y-axis.

data={{51.6,10},{16.52,18},{48.835,12},{72.975,7},{61.085,19},{78.362,11},{24.312,25},
{26.818,5},{37.336,5},{70.647,18},{48.783,35},{40.379,35},{75.926,8},{23.226,8},
{54.727,18},{116.046,12},{20.597,20},{21.973,28},{74.821,7},{56.408,51}};


After rescaling and truncating the distribution, it looks like this:

dnorm=Transpose[Rescale/@Transpose[data]];
pdf[x_,y_]=PDF[TruncatedDistribution[{{0,∞},{0,∞}},
ContourPlot[pdf[x,y],{x,0,1},{y,0,1}]


I am interested in the profile across a line parallel to the y-axis (say e.g. $x=0.5$) and do handle this profile as a distribution for this specific x-value. I need the 95% quantile for this specific distribution at the specific x-value.

Finally, I am interested in all the quantile-values for all x-values between 0 and 1 (in the rescaled domain) and do rescale them back afterwards.

I have written the following function, which works but need approx. 10 seconds (even without the line plot and sampling at 20 data points between 0 and 1). I have to evaluate hundreds of data sets with this function, so I need to speed it up.

myquantile[d_List,q_]:=Block[{dnorm,pdf,x,y,func,int,ξ,ψ,sol,invx,invy,lp},
dnorm=Transpose[Rescale/@Transpose[d]];
pdf[x_,y_]=PDF[TruncatedDistribution[{{0,∞},{0,∞}},
func[ξ_]:=Block[{},
int=NIntegrate[pdf[ξ,y],{y,0,1}];
ψ/.FindRoot[NIntegrate[pdf[ξ,y],{y,0,ψ}]/int==q,{ψ,0.8,0.9},Evaluated->False]];
{invx[x_],invy[x_]}=Rescale[x,{0,1},#]&/@({Min[#],Max[#]}&/@Transpose[d]);
sol=Interpolation[Table[{invx[x],invy[func[x]]},{x,0,1,0.05}]
];
lp=Plot[sol[x],{x,invx[0],invx[1]},AxesOrigin->{invx[0],invy[0]}];
{sol,lp}
];
myquantile[data, 0.95] // Timing


The final InterpolatingFunction looks like this (in the regular domain):

Can anybody point me in the right direction? Is it possible to get the line distribution in a different and faster way? A MarginalDistribution is pretty much what I want, but it is the distribution for all x-values and it gives me no information on individual x-values.

Thanks!

-
@0x4A4D: how have you changed the [Psi] etc. into the "real" symbols? It was showing up differently, when I posted this question? –  akm Jun 20 '13 at 13:21
Well, magic... okay, no. I knew the Unicode for these, so I put 'em in. –  Ｊ. Ｍ. Jun 20 '13 at 13:24

You can setup a partial differential equation to integrate along $y$. It costs less than 0.1 second on my old machine.

sol = NDSolve[{
D[cdf[x, y], y] == pdf[x, y],
cdf[x, 0] == 0
},
cdf, {x, 0, 1}, {y, 0, 1}]

quantileLine = ContourPlot[Evaluate[
cdf[x, y]/cdf[x, 1] == .9 /. sol[[1]]],
{x, 0, 1}, {y, 0, 1}]


Or if you want an explicit InterpolatingFunction:

quantileLineFunc = Interpolation@
Cases[quantileLine,
GraphicsComplex[pts_, others__] :>
Part[pts,
Cases[others, Line[idx_] :> idx, ∞][[1]]
]
][[1]]

-
@ Silvia: thanks a lot; your answer is very obvious, since you just use the dependency between the cdf and pdf. I was completely blocked in mind and haven't thought about other directions to solve this. Again, thanks a lot! –  akm Jun 20 '13 at 15:31
@akm You're welcome! –  Silvia Jun 20 '13 at 15:43