Log of Integral

**chogo** · July 21st 2010, 08:39 AM

I posted this problem in a different form a year ago! I have encountered it again and would value any suggestions on what approach i should take or what i should read into

This is a part of my equation. i need to evaluate this

$f(x) = \int_m^y (1-x)^N*(x)^N*e^{Nsx} dx$

where N and S are constants

my problem is that while i can use numerical integration, if the value of N is very large, then i cannot evaluate the exponent (computers can not have integers greater than 1x10^301)

so i want to take the log of that function. That would be trivial without the integral, but i am not sure how to do it with the integral.

How do i take logs? Thanks for ANY help

Chogo

**Ackbeet** · July 21st 2010, 09:11 AM

Do you mean $f(y)=...$ ? When you perform the integration, there are no x's left.

**chogo** · July 21st 2010, 09:58 AM

I mean $f(m<x<n|N,s) = \int_m^n (1-x)^N*(x)^N*e^{Nsx} dx$

**Ackbeet** · July 21st 2010, 11:32 AM

Well, if you were to look at the indefinite integral, you could differentiate it and then take the log of both sides. You'd get

$\ln(f'(y))=N\ln(1-y)+N\ln(y)+N s y.$

However, you're ultimately interested in the function f. It's possible to differentiate this equation to obtain the following:

$\displaystyle{\frac{f''(y)}{f'(y)}=-\frac{N}{1-y}+\frac{N}{y}+N s.}$

You could then solve this differential equation numerically with the proper initial conditions. Doing this might be more numerically stable than your integration.

I have no idea if this will work. It's just an idea.

**chogo** · July 21st 2010, 12:24 PM

Many thanks

What is this type of approach called? so i can read up on it

this looks an interesting approach, what happened to the integration bounds?

chogo

**Ackbeet** · July 21st 2010, 12:26 PM

I suppose it's related to logarithmic differentiation. I'm not sure this particular method has a name, though. Let me know how it works out for you.

[EDIT]: The integration bounds went away, because in order to take derivatives, I had to have an indefinite integral. You can probably recover the limits by choosing an appropriate initial condition, and probably stopping your numerical solution of the DE at the appropriate location.

**chogo** · July 23rd 2010, 11:07 AM

Firstly thanks for your help, I will be sure to acknowledge you fully in the publication. I spoke to a colleague and we tried to use your approach. Let me reformulate the question for you in full

this is the function $P(x) = \frac{\int_m^n (1-x)^{N-1} x^{N-1} e^{NSx}} {\int_0^1 (1-x)^{N-1} x^{N-1} e^{NSx}}$

I need to evaluate this integral (before i gave you just the top part for simplicity, but this is the whole function). The problem is the SAME as before, i need to take logs as i cant evaluate the exponent

so Let me rewrite this as

$P(x) = \frac{\int_m^n f(x) dx}{\int_0^1 g(x) dx}$

now using your approach of logarithmic differentiation

$=\frac{\int_m^n log(f(x)) dx}{\int_0^1 log(g(x)) dx}$

Now the key part. I have tried something different (let me know what you think) Consider the numerator

$log f(x) = \int \frac{f'}{f} dx$

and therefore

$\frac{f'}{f} = \frac{d}{dx} log(f(x)$

if i can find f and f' (which i have done - its long but not hard)

$f=(1-x)^{N-1} x^{N-1} e^{Nsx}$

and

$f' = x^{N-1} e^{NSx} (1-x)^{N-1} (2N-2+NS)$

so now by integrating f/f' i can say

$log f(x) = NSx + 2Nx - 2x + C$

and so

$P(x) =\frac{\int_m^n log(f(x)) dx}{\int_0^1 log(g(x)) dx} = \frac{\int_m^n e^{ NSx + 2Nx - 2x}}{\int_0^1 e^{ NSx + 2Nx - 2x}}$

now i am stuck. I have greatly simplified the expression, but i dont know what do do next. Anyone have any ideas? I still have the exponent

I wanted to go with your changing the integral to a differential equation route but with the full equation i was unable to do it.

This is a lot to read, and i am asking for alot of help. So I will understand if you dont reply or have the time for this. Also apologies for my shoddy mathematics.

**Ackbeet** · July 23rd 2010, 04:43 PM

Hmm. Unfortunately, you've used several steps that are not permissible mathematically. The goal is to find $P$ , which is not a function of x, since x is the variable of integration. It is a function of m and n, however.

We take each integral separately. We examine the first one, and strip away the limits of integration to examine the indefinite integral:

$\int(1-x)^{N-1} x^{N-1} e^{NSx}\,dx=f(x).$

We should note that the goal, then, of this problem is to evaluate

$P=\frac{f(n)-f(m)}{f(1)-f(0)}.$ This assumes that the function f is continuous.

Taking the derivative of both sides yields

$f'(x)=(1-x)^{N-1} x^{N-1} e^{NSx}.$ We take the logarithm of both sides to yield

$\ln(f'(x))=(N-1)\ln(1-x)+(N-1)\ln(x)+N s x.$

Finally, we take one more derivative to yield

$\displaystyle{\frac{f''(x)}{f'(x)}=\frac{N-1}{1-x}+\frac{N-1}{x}+N s}.$

Solving this numerically shouldn't be difficult in principle. The actual difficulty is in getting initial conditions. You can get the $f'(x)$ values easily enough from the equation above. But the $f(x)$ values are precisely what we are trying to find!

Here is what I recommend: a modified "shooting" algorithm. Normally, when solving a second-order ODE, you need initial conditions on $f(x)$ and $f'(x).$ You have the second, but not the first. However, you do know where you want your algorithm to end up. You can easily evaluate $f'(n), f'(m), f'(1),$ and $f'(0).$ So, what you do is you take the lower limit of the upper integral, n, as your starting point. You get your $f'(n)$ for one initial condition. You guess through a whole range of initial values for $f(n).$ Whichever one of those values gets you to m with the derivative being $f'(m)$ was the correct initial condition for $f(x),$ and is, incidentally, the value of $f(n)$ . Using that numerical solution, you can evaluate to find $f(m)$ . Repeat for the $0$ to $1$ case. You see how this works?

I think this algorithm might just work for you.

One more comment: I would double-check your exponents. They seem to have changed from the original post.

**chogo** · July 24th 2010, 06:18 AM

Thank you. I will try to solve this in MATLAB

i will let you know how it works. Thanks again for your time

**Ackbeet** · July 24th 2010, 07:17 AM

Just as a warning: the 0 to 1 case might give you some convergence trouble. On the other hand, perhaps the original integral there isn't so unstable? If I look at

$\int_0^1 (1-x)^{N-1} x^{N-1} e^{NSx},$

because the values of x and 1-x are between 0 and 1, I think your exponentiations are going to be smaller.

What are the nominal sizes for n, m, N, and S? What interval of values do they range over?

**chogo** · July 24th 2010, 07:34 AM

Hey Ackbeet

I was just trying to progam this up now. Is there a more conventional way of writing the second order DE, ive never seen one in this form and am trying to use the ODE45 solver in MATLAB

0 and 1 will definitely give convergence trouble. The integral is undefined for limit $0$ if $N-1<0$ as the term $x^{N-1} -> 0^{N-1}$ which is undefined

But I am modelling a biological system and so i can limits 0.001 and 0.999 instead of 0 and 1. This should help i think

these are the interval ranges

N - 5-10,000
s=0.7-2.4
n,m are bin sizes for the integral typically 0.01 or 0.001 with the range {0.001, 0.999} - the same as above

Thanks Again

Samir Bhatt

**Ackbeet** · July 24th 2010, 08:44 AM

You can write the DE as

$\displaystyle{f''(x)=f'(x)\left[\frac{N-1}{1-x}+\frac{N-1}{x}+N s\right]},$ or

$\displaystyle{f''(x)-f'(x)\left[\frac{N-1}{1-x}+\frac{N-1}{x}+N s\right]=0}.$

The latter is very conventional.

I would concur with your approach to handling 0 and 1. Just go a little inside those boundaries, and hopefully you'll be ok.

Incidentally, I gave the following integral to Mathematica, increased the recursion depth to 15, and got the following result with no complaints:

$\int_{0}^{1}(1-x)^{10000-1}x^{10000-1}e^{10000*2.4*x}\,dx,$

with the Mathematica command:

NIntegrate[(1 - x)^(10000 - 1)*x^(10000 - 1)*Exp[10000*2.4*x], {x, 0, 1},
MaxRecursion -> 15],

yielded

$6.5512155712\times 10^{554}.$

If that's the worst-case scenario, then maybe you just need to switch over to Mathematica!

I should also mention that Mathematica could easily compute

$e^{10000*2.4}=1.16833040227\times 10^{10423},$

considerably larger than the limit you mentioned.

I'm running Mathematica Version 4 on Windows XP Professional, in case you want to know.

Thread: Log of Integral

LinkBack

Thread Tools

Display

Log of Integral

Similar Math Help Forum Discussions

Use Cauchy's integral theorem for derivatives to evaluate the contour integral.

Graph the region of integration and evaluate the double integral (integral from 0 to

Explaining result: comparing line integral and two-form integral

[SOLVED] Line integral, Cauchy's integral formula

write [integral of] dx/SQRT(sinx) in terms of an elliptic integral

Search Tags