[R] Computational problems in R

Steven McKinney smckinney at bccrc.ca
Tue Oct 28 20:06:08 CET 2008


I misinterpreted the suggestion from Duncan
Murdoch, and I think Amy did too, judging
from her comment about B possibly being negative.

Duncan Murdoch kindly explained the idea
off-list, many thanks to him for taking the time.
To avoid further confusion, I've set up a basic
example explaining the suggestion.

Suppose this is the sum 
sum(exp(c_i/d)) i = 1, ..., n
that needs evaluating:

exp(-71) + exp(12) + exp(842) + exp(2100) + exp(141) + exp(-710)

The maximum term in the sum is exp(2100) 
so in the notation of the identity
A+B = A*(1 + B/A) = exp(log(A) + log(1 + B/A))

A = exp(2100)
B = exp(-71) + exp(12) + exp(842) + exp(141) + exp(-710)

so sum(exp(c_i/d)) i = 1, ..., n is
exp(2100) * ( 1 + ((exp(-71) + exp(12) + exp(842) + exp(141) + exp(-710))/exp(2100)) )
= exp(2100) * ( 1 + exp(-71 - 2100) + exp(12 - 2100) + exp(842 - 2100) + exp(141 - 2100) + exp(-710 - 2100) )

and of course all the terms 
exp(-71 - 2100) 
exp(12 - 2100) 
exp(842 - 2100) 
exp(141 - 2100) 
exp(-710 - 2100)
are very close to zero, so the sum
sum(exp(c_i/d)) i = 1, ..., n is 
exp(2100)
or 2100 on the log scale.

The terms in B are of course never negative as they
are all exponential quantities greater than zero,
so that's not an issue in this calculation (this
'B' not being the same as the 'B' in Amy's original
description).

Thanks again to Duncan Murdoch, and my apologies
for my misinterpretation.


Steven McKinney


-----Original Message-----
From: Robin Hankin [mailto:rksh1 at cam.ac.uk]
Sent: Tue 10/28/2008 12:58 AM
To: A.Noufaily
Cc: Duncan Murdoch; r-help at r-project.org; markleeds at verizon.net; Steven McKinney; Xiaoxu LI
Subject: Re: [R] Computational problems in R
 
Hello.

The Brobdingnag package uses that identity for a logarithmic
representation and also has a hack for negative numbers.


HTH

rksh



A.Noufaily wrote:
>  
> Many thanks for your suggestions...
>
> I am still checking which one is the most useful for my simulations.
>
> Concerning using logs, this might be very helpful, but I am not sure if
> I can use the following:
> A+B = A*(1 + B/A)
>     = exp(log(A) + log(1 + B/A))
> because unfortunately B can be negative.
> However, I might still use logs in case (1 + B/A)>0.
>
> Regards,
>  
> Amy
>
> -----Original Message-----
> From: Duncan Murdoch [mailto:murdoch at stats.uwo.ca] 
> Sent: Saturday, October 25, 2008 11:36 AM
> To: Steven McKinney
> Cc: A.Noufaily; r-help at r-project.org
> Subject: Re: [R] Computational problems in R
>
> On 24/10/2008 9:50 PM, Steven McKinney wrote:
>   
>> I suspect there's a deeper issue here.
>> sum(exp(yi)) when large yi occur is
>> problematic.  exp(yi) for yi>710 is
>> just a huge number, and summing additional values only makes the 
>> overall sum larger as all components of the summation are positive.
>> There's no way around that.
>>     
>
> Sure there is, and you quoted it below.  Work on a log scale.  The log
> of exp(yi) is yi, and it sounds as though the yi values are manageable.
>
> You might end up knowing that the log of the final answer is 20000 and
> not be able to evaluate exp(20000) in R, but you still know that the
> answer is exp(20000).
>
> Duncan Murdoch
>   
>> You could try this with Robin Hankins'
>> package "brobdingnag" which can handle bunches of bizarrely large 
>> numbers.
>>
>> What kind of process are you studying?
>> What kind of process generates values
>> such as exp(8/0.01) when other values
>> are much smaller?  Are these outliers
>> in an otherwise well-behaved
>> data set?  Perhaps then they need
>> to be set aside and investigated
>> separately, etc.
>>
>>
>> Steven McKinney
>>
>> Statistician
>> Molecular Oncology and Breast Cancer Program British Columbia Cancer 
>> Research Centre
>>
>> email: smckinney +at+ bccrc +dot+ ca
>>
>> tel: 604-675-8000 x7561
>>
>> BCCRC
>> Molecular Oncology
>> 675 West 10th Ave, Floor 4
>> Vancouver B.C. 
>> V5Z 1L3
>> Canada
>>
>>
>>
>>
>> -----Original Message-----
>> From: r-help-bounces at r-project.org on behalf of Duncan Murdoch
>> Sent: Fri 10/24/2008 4:04 PM
>> To: A.Noufaily
>> Cc: r-help at r-project.org
>> Subject: Re: [R] Computational problems in R
>>  
>> On 24/10/2008 12:42 PM, A.Noufaily wrote:
>>     
>>> Dear all,
>>>
>>> I would be grateful if anyone can help me with the following:
>>>
>>> My aim is to compute explicitely the sum S=A+B where 
>>> A=sum(exp(c_i/d)), i=1,...,n; B, c_i, and d are real numbers with 
>>> -Inf<B,c_i<+Inf; and d>0.
>>> The problem is that when c_i/d >710 (for some i) R is setting
>>> exp(c_i/d) to be equal to +Inf and hence the whole summation S.
>>> So in simple cases where for example c_i=8 (for some i), and d=0.01 
>>> the whole summation is turning out to be infinite.
>>> Is there a way to get round that in R?
>>> Can anyone suggest any computational trick to calculate S when 
>>> c_i/d>710 (for some i)?
>>>       
>> Work on a log scale.  Use the identity that
>>
>> A+B = A*(1 + B/A)
>>      = exp(log(A) + log(1 + B/A))
>>
>> (where you chose A to be the biggest term in the sum).
>>
>> Duncan Murdoch
>>
>>     
>>> Any suggestions would be much appreciated.
>>>
>>> Regards,
>>>
>>> Amy
>>>
>>>
>>>
>>>
>>>
>>> ---------------------------------
>>> The Open University is incorporated by Royal Charter (RC 000391), an
>>>       
> exempt charity in England & Wales and a charity registered in Scotland
> (SC 038302).
>   
>>> 	[[alternative HTML version deleted]]
>>>
>>> ______________________________________________
>>> R-help at r-project.org 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.
>>>       
>> ______________________________________________
>> R-help at r-project.org 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.
>>
>> ______________________________________________
>> R-help at r-project.org 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.
>>     
>
>
> ---------------------------------
> The Open University is incorporated by Royal Charter (RC 000391), an exempt charity in England & Wales and a charity registered in Scotland (SC 038302).
>
> ______________________________________________
> R-help at r-project.org 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.
>   



More information about the R-help mailing list