[Rd] Computing means, variances and sums

2006-02-19 Thread Prof Brian Ripley
There has been a recent thread on R-help on this, which resurrected 
concepts from bug reports PR#1228 and PR#6743.  Since the discussion has 
included a lot of erroneous 'information' based on misunderstandings of 
floating-point computations, this is an attempt to set the record straight 
and explain the solutions adopted.

The problem was that var(rep(0.02, 10)) was observed to be (on some 
machines) about 1e-35.  This can easily be explained.

0.02 is not an exactly repesented binary fraction.  Repeatedly adding the 
represented quantity makes increasing rounding errors relative to the 
exact computation, so (on Sparc Solaris)

> var(rep(0.02, 10))
[1] 1.337451e-35
> options(digits=18)
> sum(rep(0.02, 10))
[1] 0.19998
> sum(rep(0.02, 10)) -0.2
[1] -2.775557561562891e-17
> sum(rep(0.02, 10))/10 -0.02
[1] -3.469446951953614e-18
> z <- sum(rep(0.02, 10))/10 -0.02
> 10*z^2/9
[1] 1.3374513502689138e-35

and so the non-zero variance is arising from (x[i] - mean) being non-zero.
(I did check that was what was happening at C level.)

There has been talk of other ways to arrange these computations, for 
example Kahan summation and Welford's algorithm (see Chan & Lewis, 1979, 
CACM 22, 526-531 and references therein).  However, R already used the 
two-pass algorithm which is the most accurate (in terms of error bounds) 
in that reference.

Why are most people seeing 0?  Because the way computation is done in 
modern FPUs is not the computation analysed in early numerical analysis 
papers, including in Chan & Lewis.  First, all FPUs that I am aware of 
allow the use of guard digits, effectively doing intermediate computations 
to one more bit than required.  And many use extended precision registers 
for computations which they can keep in FP registers, thereby keeping 
another 10 or more bits.  (This includes R on most OSes on ix86 CPUs, the 
exception being on Windows where the FPU has been reset by some other 
software. Typically it is not the case for RISC CPUs, e.g. Sparc.)

The use of extended precision registers invalidates the textbook 
comparisons of accuracy in at least two ways.  First, the error analysis 
is different if all intermediate results are stored in extended precision. 
Second, the simpler the algorithm, the more intermediate results which 
will remain in extended precision.  This means that (for example) Kahan 
summation is usually less accurate than direct summation on real-world 
FPUs.

There are at least two steps which can be done to improve accuracy.
One is to ensure that intermediate results are stored in extended 
precision.  Every R platform of which I am aware has a 'long double' type 
which can be used.  On ix86 this is the extended precision type used 
internally in the FPU and so the cost is zero or close to zero, whereas on 
a Sparc the extra precision is more but there is some cost.  The second 
step is to use iterative refinement, so the final part of mean.default 
currently is

 ## sum(int) can overflow, so convert here.
 if(is.integer(x)) x <- as.numeric(x)
 ## use one round of iterative refinement
 res <- sum(x)/n
 if(is.finite(res)) res + sum(x-res)/n else res

This is a well-known technique in numerical linear algebra, and improves 
the accuracy whilst doubling the cost.  (This is about to be replaced by 
an internal function to allow the intermediate result to be stored in a 
long double.)

Note the is.finite(res) there.  R works with the extended IEC60059 (aka 
IEEE754) quantities of Inf, -Inf and NaN (NA being a special NaN).  The 
rearranged computations do not work correctly for those quantities.  So 
although they can be more 'efficient' (in terms of flops), they have to be 
supplemented by some other calculation to ensure that the specials are 
handled correctly.  Remember once again that we get both speed and 
accuracy advantages by keeping computations in FPU registers, so 
complicating the code has considerable cost.

R-devel will shortly use long doubles for critical intermediate results 
and iterative refinement for calculations of means.  This may be slower 
but it would be an extreme case in which the speed difference was 
detectable.  Higher accuracy has a cost too: there are several packages 
(dprep and mclust are two) whose results are critically dependent on fine 
details of computations and will for example infinite-loop if an optimized 
BLAS is used on AMD64.


The choice of algorithms in R is an eclectic mixture of accuracy and 
speed.  When (some of) R-core decided to make use of a BLAS for e.g. 
matrix products this produced a large speed increase for those with 
optimized BLAS (and a small speed decrease for others), but it did result 
in lower accuracy and problems with NAs etc (and the alternative 
algorithms have since been added back to cover such cases).  But it seems 
that nowadays few R users understand the notion of rounding error, and it 
is easier to make the computations maximally accurate than to keep 
e

Re: [Rd] Bug in Sweave? -- scoping problem? (PR#8615)

2006-02-19 Thread berwin
G'day Duncan,

> "DM" == murdoch  <[EMAIL PROTECTED]> writes:

DM> I have found a strange scoping problem in Sweave.  [...]

DM> The strange thing is that while the value in partytotal is
DM> output correctly as [...]

DM> but the dotchart contains the wrong values: it shows sorted
DM> values, but not sorted names, [...]
No bug, but a feature. :)
(But then, a well known software producer seems to use `feature' as a
euphemism for `bug', so perhaps I shouldn't call it a feature.)

There was recently a discussion on r-help (?) about Sweave producing
different output in the text and plot when random numbers were
generated, and it seems as if you have run into the same trap:

  Code in chunks that produce pictures is executed several times.
  First, to produce the output in the text.  And then once more
  for *each* format in which the figure has to be produced.  I.e.,
  if you want a PDF and a PostScript version of the figure, the
  code is executed a total of three times.

All instances of this feature reported so far involved commands that
produced random numbers and the poster was surprised that the output
in the text and the figures differed (and that the two figures were
different).

In your case the first execution of the code assigns the names to
partytotal, sorts partytotal and produces the data.  On the next
execution, when the picture is produced, partytotal is already sorted
but you reassign the names.  Then the partytotal is sorted again and
the plot is produced.  But that re-assigning of names lead to the
disconnect between values and names.

Hope this helps.

Cheers,

Berwin

== Full address 
Berwin A Turlach  Tel.: +61 (8) 6488 3338 (secr)   
School of Mathematics and Statistics+61 (8) 6488 3383 (self)  
The University of Western Australia   FAX : +61 (8) 6488 1028
35 Stirling Highway   
Crawley WA 6009e-mail: [EMAIL PROTECTED]
Australiahttp://www.maths.uwa.edu.au/~berwin

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Bug in Sweave? -- scoping problem? (PR#8615)

2006-02-19 Thread murdoch
Berwin A Turlach wrote:
> G'day Duncan,
> 
> 
>>"DM" == murdoch  <[EMAIL PROTECTED]> writes:
> 
> 
> DM> I have found a strange scoping problem in Sweave.  [...]
> 
> DM> The strange thing is that while the value in partytotal is
> DM> output correctly as [...]
> 
> DM> but the dotchart contains the wrong values: it shows sorted
> DM> values, but not sorted names, [...]
> No bug, but a feature. :)
> (But then, a well known software producer seems to use `feature' as a
> euphemism for `bug', so perhaps I shouldn't call it a feature.)
> 
> There was recently a discussion on r-help (?) about Sweave producing
> different output in the text and plot when random numbers were
> generated, and it seems as if you have run into the same trap:
> 
>   Code in chunks that produce pictures is executed several times.
>   First, to produce the output in the text.  And then once more
>   for *each* format in which the figure has to be produced.  I.e.,
>   if you want a PDF and a PostScript version of the figure, the
>   code is executed a total of three times.

Thanks, that's what caught me.  Is that a quote from the discussion, or 
from the docs somewhere?  It makes sense in hindsight, but it's not 
obvious ahead of time, so it should be stated fairly prominently in the 
docs.

Duncan Murdoch
> 
> All instances of this feature reported so far involved commands that
> produced random numbers and the poster was surprised that the output
> in the text and the figures differed (and that the two figures were
> different).
> 
> In your case the first execution of the code assigns the names to
> partytotal, sorts partytotal and produces the data.  On the next
> execution, when the picture is produced, partytotal is already sorted
> but you reassign the names.  Then the partytotal is sorted again and
> the plot is produced.  But that re-assigning of names lead to the
> disconnect between values and names.
> 
> Hope this helps.
> 
> Cheers,
> 
> Berwin
> 
> == Full address 
> Berwin A Turlach  Tel.: +61 (8) 6488 3338 (secr)   
> School of Mathematics and Statistics+61 (8) 6488 3383 (self)  
> The University of Western Australia   FAX : +61 (8) 6488 1028
> 35 Stirling Highway   
> Crawley WA 6009e-mail: [EMAIL PROTECTED]
> Australiahttp://www.maths.uwa.edu.au/~berwin

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Error message while installing quatreg in ox s (PR#8616)

2006-02-19 Thread stefano iacus
This is not an R bug, nor it realtes to quantreg. It is not even an  
error but a warning.
It is due to the fact that you don't have enough writing permission  
and it is related the R.app way of handling package installation.
stefano

Il giorno 18/feb/06, alle ore 21:07, [EMAIL PROTECTED] ha  
scritto:

> Full_Name: Alok Krishen
> Version: 2.2.1
> OS: OS X
> Submission from: (NULL) (68.221.92.169)
>
>
> When install.packages("quantreg") produces the following error message
> cannot create HTML package index in: make.packages.html()
>
> __
> R-devel@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Computing means, variances and sums

2006-02-19 Thread Spencer Graves
Dear Professors Ripley & Murdoch, & others:

  If this were just an issue of computations like var(rep(0.02, 10)) 
producing something other than 0 on certain platforms (e.g., 
combinations of operating systems and microprocessors), then I would 
suggest it be documented as an FAQ and left as a tool to help explain 
finite precision arithmetic and rounding issues to people who don't yet 
understand those concepts.

  However, Duncan's comment shows that it is more than that, namely 
that certain DLLs change the precision of the fpu (floating point unit?) 
from 64 to 53 bit matissas AND DON'T RESET THEM.  When this is not 
detected, it could make the difference between a useable answer and 
nonsense in poorly conditioned applications where those 9 bits might be 
important.  For me, that problem is NOT that one occasionally gets 
nonsense from a poorly conditioned compution.  Rather it is that the 
SAME computation could give a useful answer in one case and nonsense a 
few seconds later ON THE SAME COMPUTER, operating system, etc.

  To test my understanding, I simplified Barry Zajdik's example further:

 > var(rep(.2, 3))
[1] 0
 > RSiteSearch("convert to binary")
A search query has been submitted to http://search.r-project.org
The results page should open in your browser shortly
 > var(rep(.2, 3))
[1] 1.18e-33
 > sessionInfo()
R version 2.2.1, 2005-12-20, i386-pc-mingw32

attached base packages:
[1] "methods"   "stats" "graphics"  "grDevices" "utils" "datasets"
[7] "base"

  This indicates there is a problem that perhaps should eventually be 
fixed.  However, please do NOT spend time on this because I suggested it 
was an issue.  The conditions under which this would create problems for 
anyone are still so rare that I would not want to alter anyone's work 
priorities for it.

  spencer graves
p.s.  If my computations are correct, 0.2 = 0*/2 + 0/4 + 1/8 + 1/16 + 
0/32 + 0/64 + 1/128 + 1/256 + 0/512 + 0/1024 + 1/2048 + 1/4096 + ... = 
0.3h.  Perhaps someone can extend this to an FAQ to help 
explain finite precision arithmetic and rounding issues.

Prof Brian Ripley wrote:
> There has been a recent thread on R-help on this, which resurrected 
> concepts from bug reports PR#1228 and PR#6743.  Since the discussion has 
> included a lot of erroneous 'information' based on misunderstandings of 
> floating-point computations, this is an attempt to set the record 
> straight and explain the solutions adopted.
> 
> The problem was that var(rep(0.02, 10)) was observed to be (on some 
> machines) about 1e-35.  This can easily be explained.
> 
> 0.02 is not an exactly repesented binary fraction.  Repeatedly adding 
> the represented quantity makes increasing rounding errors relative to 
> the exact computation, so (on Sparc Solaris)
> 
>> var(rep(0.02, 10))
> 
> [1] 1.337451e-35
> 
>> options(digits=18)
>> sum(rep(0.02, 10))
> 
> [1] 0.19998
> 
>> sum(rep(0.02, 10)) -0.2
> 
> [1] -2.775557561562891e-17
> 
>> sum(rep(0.02, 10))/10 -0.02
> 
> [1] -3.469446951953614e-18
> 
>> z <- sum(rep(0.02, 10))/10 -0.02
>> 10*z^2/9
> 
> [1] 1.3374513502689138e-35
> 
> and so the non-zero variance is arising from (x[i] - mean) being non-zero.
> (I did check that was what was happening at C level.)
> 
> There has been talk of other ways to arrange these computations, for 
> example Kahan summation and Welford's algorithm (see Chan & Lewis, 1979, 
> CACM 22, 526-531 and references therein).  However, R already used the 
> two-pass algorithm which is the most accurate (in terms of error bounds) 
> in that reference.
> 
> Why are most people seeing 0?  Because the way computation is done in 
> modern FPUs is not the computation analysed in early numerical analysis 
> papers, including in Chan & Lewis.  First, all FPUs that I am aware of 
> allow the use of guard digits, effectively doing intermediate 
> computations to one more bit than required.  And many use extended 
> precision registers for computations which they can keep in FP 
> registers, thereby keeping another 10 or more bits.  (This includes R on 
> most OSes on ix86 CPUs, the exception being on Windows where the FPU has 
> been reset by some other software. Typically it is not the case for RISC 
> CPUs, e.g. Sparc.)
> 
> The use of extended precision registers invalidates the textbook 
> comparisons of accuracy in at least two ways.  First, the error analysis 
> is different if all intermediate results are stored in extended 
> precision. Second, the simpler the algorithm, the more intermediate 
> results which will remain in extended precision.  This means that (for 
> example) Kahan summation is usually less accurate than direct summation 
> on real-world FPUs.
> 
> There are at least two steps which can be done to improve accuracy.
> One is to ensure that intermediate results are stored in extended 
> precision.  Every R platform of which I am aware has a 'long double' 
> type which can be

Re: [Rd] Computing means, variances and sums

2006-02-19 Thread hadley wickham
> p.s.  If my computations are correct, 0.2 = 0*/2 + 0/4 + 1/8 + 1/16 +
> 0/32 + 0/64 + 1/128 + 1/256 + 0/512 + 0/1024 + 1/2048 + 1/4096 + ... =
> 0.3h.  Perhaps someone can extend this to an FAQ to help
> explain finite precision arithmetic and rounding issues.

This is drifting a bit off topic, but the other day I discovered this
rather nice illustration of the perils of finite precision arithmetic
while creating a contrast matrix:

> n <- 13
> a <- matrix(-1/n, ncol=n, nrow=n) + diag(n)
> rowSums(a)
 [1]  2.775558e-16  2.775558e-16  5.551115e-17  5.551115e-17  5.551115e-17
 [6]  5.551115e-17  0.00e+00 -5.551115e-17  0.00e+00  5.551115e-17
[11]  1.110223e-16  1.665335e-16  2.220446e-16

Not only do most of the rows not sum to 0, they do not even sum to the
same number!  It is hard to remember the familiar rules of arithmetic
do not always apply.

Hadley

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Computing means, variances and sums

2006-02-19 Thread Prof Brian Ripley
As I have said before on R-help, it is _already_ fixed in this (and 
related cases) and the reasons it is fixed were explained here.  If you 
read carefully, you will have noticed that some platforms (such as Sparc) 
only use 53 bits (+ a guard bit).

Not that I accept that the 53-bit calculation is `nonsense', or anything 
close to it.

On Sun, 19 Feb 2006, Spencer Graves wrote:

> Dear Professors Ripley & Murdoch, & others:
>
> If this were just an issue of computations like var(rep(0.02, 10)) 
> producing something other than 0 on certain platforms (e.g., combinations of 
> operating systems and microprocessors), then I would suggest it be documented 
> as an FAQ and left as a tool to help explain finite precision arithmetic and 
> rounding issues to people who don't yet understand those concepts.
>
> However, Duncan's comment shows that it is more than that, namely 
> that certain DLLs change the precision of the fpu (floating point unit?) from 
> 64 to 53 bit matissas AND DON'T RESET THEM.  When this is not detected, it 
> could make the difference between a useable answer and nonsense in poorly 
> conditioned applications where those 9 bits might be important.  For me, that 
> problem is NOT that one occasionally gets nonsense from a poorly conditioned 
> compution.  Rather it is that the SAME computation could give a useful answer 
> in one case and nonsense a few seconds later ON THE SAME COMPUTER, operating 
> system, etc.
>
> To test my understanding, I simplified Barry Zajdik's example 
> further:
>
>> var(rep(.2, 3))
> [1] 0
>> RSiteSearch("convert to binary")
> A search query has been submitted to http://search.r-project.org
> The results page should open in your browser shortly
>> var(rep(.2, 3))
> [1] 1.18e-33
>> sessionInfo()
> R version 2.2.1, 2005-12-20, i386-pc-mingw32
>
> attached base packages:
> [1] "methods"   "stats" "graphics"  "grDevices" "utils" "datasets"
> [7] "base"
>
> This indicates there is a problem that perhaps should eventually be 
> fixed.  However, please do NOT spend time on this because I suggested it was 
> an issue.  The conditions under which this would create problems for anyone 
> are still so rare that I would not want to alter anyone's work priorities for 
> it.
>
> spencer graves
> p.s.  If my computations are correct, 0.2 = 0*/2 + 0/4 + 1/8 + 1/16 + 0/32 + 
> 0/64 + 1/128 + 1/256 + 0/512 + 0/1024 + 1/2048 + 1/4096 + ... = 
> 0.3h.  Perhaps someone can extend this to an FAQ to help explain 
> finite precision arithmetic and rounding issues.
>   Prof Brian Ripley wrote:
>> There has been a recent thread on R-help on this, which resurrected 
>> concepts from bug reports PR#1228 and PR#6743.  Since the discussion has 
>> included a lot of erroneous 'information' based on misunderstandings of 
>> floating-point computations, this is an attempt to set the record straight 
>> and explain the solutions adopted.
>> 
>> The problem was that var(rep(0.02, 10)) was observed to be (on some 
>> machines) about 1e-35.  This can easily be explained.
>> 
>> 0.02 is not an exactly repesented binary fraction.  Repeatedly adding the 
>> represented quantity makes increasing rounding errors relative to the exact 
>> computation, so (on Sparc Solaris)
>> 
>>> var(rep(0.02, 10))
>> 
>> [1] 1.337451e-35
>> 
>>> options(digits=18)
>>> sum(rep(0.02, 10))
>> 
>> [1] 0.19998
>> 
>>> sum(rep(0.02, 10)) -0.2
>> 
>> [1] -2.775557561562891e-17
>> 
>>> sum(rep(0.02, 10))/10 -0.02
>> 
>> [1] -3.469446951953614e-18
>> 
>>> z <- sum(rep(0.02, 10))/10 -0.02
>>> 10*z^2/9
>> 
>> [1] 1.3374513502689138e-35
>> 
>> and so the non-zero variance is arising from (x[i] - mean) being non-zero.
>> (I did check that was what was happening at C level.)
>> 
>> There has been talk of other ways to arrange these computations, for 
>> example Kahan summation and Welford's algorithm (see Chan & Lewis, 1979, 
>> CACM 22, 526-531 and references therein).  However, R already used the 
>> two-pass algorithm which is the most accurate (in terms of error bounds) in 
>> that reference.
>> 
>> Why are most people seeing 0?  Because the way computation is done in 
>> modern FPUs is not the computation analysed in early numerical analysis 
>> papers, including in Chan & Lewis.  First, all FPUs that I am aware of 
>> allow the use of guard digits, effectively doing intermediate computations 
>> to one more bit than required.  And many use extended precision registers 
>> for computations which they can keep in FP registers, thereby keeping 
>> another 10 or more bits.  (This includes R on most OSes on ix86 CPUs, the 
>> exception being on Windows where the FPU has been reset by some other 
>> software. Typically it is not the case for RISC CPUs, e.g. Sparc.)
>> 
>> The use of extended precision registers invalidates the textbook 
>> comparisons of accuracy in at least two ways.  First, the error analysis is 
>> different if all intermediate res

Re: [Rd] Computing means, variances and sums

2006-02-19 Thread Prof Brian Ripley
On Sun, 19 Feb 2006, hadley wickham wrote:

>> p.s.  If my computations are correct, 0.2 = 0*/2 + 0/4 + 1/8 + 1/16 +
>> 0/32 + 0/64 + 1/128 + 1/256 + 0/512 + 0/1024 + 1/2048 + 1/4096 + ... =
>> 0.3h.  Perhaps someone can extend this to an FAQ to help
>> explain finite precision arithmetic and rounding issues.
>
> This is drifting a bit off topic, but the other day I discovered this
> rather nice illustration of the perils of finite precision arithmetic
> while creating a contrast matrix:
>
>> n <- 13
>> a <- matrix(-1/n, ncol=n, nrow=n) + diag(n)
>> rowSums(a)
> [1]  2.775558e-16  2.775558e-16  5.551115e-17  5.551115e-17  5.551115e-17
> [6]  5.551115e-17  0.00e+00 -5.551115e-17  0.00e+00  5.551115e-17
> [11]  1.110223e-16  1.665335e-16  2.220446e-16
>
> Not only do most of the rows not sum to 0, they do not even sum to the
> same number!  It is hard to remember the familiar rules of arithmetic
> do not always apply.

I think you will find this example does give all 0's in R-devel, even 
on platforms like Sparc.  But users do need to remember that computer 
arithmetic is inexact except in rather narrowly delimited cases.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Computing means, variances and sums

2006-02-19 Thread Duncan Murdoch
On 2/19/2006 3:18 PM, Prof Brian Ripley wrote:
> On Sun, 19 Feb 2006, hadley wickham wrote:
> 
>>> p.s.  If my computations are correct, 0.2 = 0*/2 + 0/4 + 1/8 + 1/16 +
>>> 0/32 + 0/64 + 1/128 + 1/256 + 0/512 + 0/1024 + 1/2048 + 1/4096 + ... =
>>> 0.3h.  Perhaps someone can extend this to an FAQ to help
>>> explain finite precision arithmetic and rounding issues.
>> This is drifting a bit off topic, but the other day I discovered this
>> rather nice illustration of the perils of finite precision arithmetic
>> while creating a contrast matrix:
>>
>>> n <- 13
>>> a <- matrix(-1/n, ncol=n, nrow=n) + diag(n)
>>> rowSums(a)
>> [1]  2.775558e-16  2.775558e-16  5.551115e-17  5.551115e-17  5.551115e-17
>> [6]  5.551115e-17  0.00e+00 -5.551115e-17  0.00e+00  5.551115e-17
>> [11]  1.110223e-16  1.665335e-16  2.220446e-16
>>
>> Not only do most of the rows not sum to 0, they do not even sum to the
>> same number!  It is hard to remember the familiar rules of arithmetic
>> do not always apply.
> 
> I think you will find this example does give all 0's in R-devel, even 
> on platforms like Sparc.  

Only until the fpu precision gets changed:

 > n <- 13
 > a <- matrix(-1/n, ncol=n, nrow=n) + diag(n)
 > rowSums(a)
  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0
 > RSiteSearch('junk')
A search query has been submitted to http://search.r-project.org
The results page should open in your browser shortly

 > n <- 13
 > a <- matrix(-1/n, ncol=n, nrow=n) + diag(n)
 > rowSums(a)
  [1]  2.775558e-16  2.775558e-16  5.551115e-17  5.551115e-17  5.551115e-17
  [6]  5.551115e-17  0.00e+00 -5.551115e-17  0.00e+00  5.551115e-17
[11]  1.110223e-16  1.665335e-16  2.220446e-16

We still need to protect against these changes.  I'll put something 
together, unless you're already working on it.

The approach I'm thinking of is to define a macro to be called in risky 
situations.  On platforms where this isn't an issue, the macro would be 
null; on Windows, it would reset the fpu to full precision.

For example, RSiteSearch causes damage in the ShellExecute call in 
do_shellexec called from browseURL, so I'd add protection there.  I 
think we should also add detection code somewhere in the evaluation loop 
to help in diagnosing these problems.

> But users do need to remember that computer 
> arithmetic is inexact except in rather narrowly delimited cases.

Yes, that too.

Duncan Murdoch

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] invalid graphics state using dev.print (fwd)

2006-02-19 Thread Paul Murrell
Hi


Martin Maechler wrote:
>>"Paul" == Paul Roebuck <[EMAIL PROTECTED]>
>>on Wed, 8 Feb 2006 15:33:11 -0600 (CST) writes:
> 
> 
> Paul> On Mon, 6 Feb 2006 18:12, Simon Urbanek wrote:
> >> On Feb 6, 2006, at 5:24 PM, Paul Roebuck wrote:
> >> 
> >>> Tried on R-Sig-Mac with no responses, but I need some kind
> >>> of answer.
> >>> [...]
> >>> Does the following work on your system?
> >> 
> >> Interesting, no, it doesn't either. For png and pdf I use
> >> Quartz + quartz.save (it produces much nicer results) so
> >> I didn't really notice, but you're right. First I thought
> >> those graphics state issues are specific to the Quartz
> >> device, but you have proven that it's not. It's in fact
> >> not even Mac-specific - I have just reproduced it on a
> >> Linux box - that's why I'm moving this to R-devel.
> 
> Paul> It's been several workdays now with no responses. Could
> Paul> someone try the last three lines of code and see if they
> Paul> get the following error message?
> 
> >> x11()
> >> plot(rnorm(10))
> >> dev.print(png)
> 
> Paul> Error in dev.copy(device = function (filename = "Rplot%03d.png", 
> width =
> Paul> 480,  :
> Paul> invalid graphics state
> 
> >> traceback()
> Paul> 6: dev.copy(device = function (filename = "Rplot%03d.png", width = 
> 480,
> Paul> height = 480, pointsize = 12, gamma = 1, colortype =
> Paul> getOption("X11colortype"),
> Paul> maxcubesize = 256, bg = "white", fonts = getOption("X11fonts"),
> Paul> res = NA)
> Paul> .Internal(X11(paste("png::", filename, sep = ""), width, height,
> Paul> pointsize, gamma, colortype, maxcubesize, bg, bg, fonts,
> Paul> res)), width = 6.98715785526809, height = 6.99452568428947)
> Paul> 5: eval(expr, envir, enclos)
> Paul> 4: eval(expr, p)
> Paul> 3: eval.parent(oc)
> Paul> 2: dev.off(eval.parent(oc))
> Paul> 1: dev.print(png)
> 
> Paul> I noticed it on OS X, and Simon on Linux. 
> 
> Yes, I can confim getting the same.
> Just on Linux though (as Simon)
> 
> I'd say this should make a ``nice little''  bug.report() 
> 
> Interestingly, replacing
> 
> dev.print(png)
> 
> by  dev.copy(png) ; dev.off()  
> 
> which is about equivalent,  *does* work and so is a workaround
> to your problem.


I think the problem is that the width and height of the PNG device is 
being taken (without regard for units) from the X11 device.  So 
approximately 7 inches square screen window gets drawn into 
approximately 7 *pixel* square PNG file and (understandably) R complains 
that there is not enough room for the plot.  Another workaround is 
something like ...

dev.print(png, width=480, height=480)

... and a fix requires making dev.print() smarter so that it figures out 
that it needs to convert width/height from inches to pixels.

Paul


> Paul> Other platforms?  WFM?
> 
> Paul> TIA
> 
> >> version
> Paul> _
> Paul> platform powerpc-apple-darwin7.9.0
> Paul> arch powerpc
> Paul> os   darwin7.9.0
> Paul> system   powerpc, darwin7.9.0
> Paul> status   Patched
> Paul> major2
> Paul> minor2.1
> Paul> year 2006
> Paul> month02
> Paul> day  01
> Paul> svn rev  37245
> Paul> language R
> 
> Paul> --
> Paul> SIGSIG -- signature too long (core dumped)
> 
> __
> R-devel@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel

-- 
Dr Paul Murrell
Department of Statistics
The University of Auckland
Private Bag 92019
Auckland
New Zealand
64 9 3737599 x85392
[EMAIL PROTECTED]
http://www.stat.auckland.ac.nz/~paul/

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] invalid graphics state using dev.print (fwd)

2006-02-19 Thread Prof Brian Ripley
On Mon, 20 Feb 2006, Paul Murrell wrote:

[...]

>>>> x11()
>>>> plot(rnorm(10))
>>>> dev.print(png)

>> Paul> Error in dev.copy(device = function (filename = "Rplot%03d.png", 
>> width =
>> Paul> 480,  :
>> Paul> invalid graphics state

> I think the problem is that the width and height of the PNG device is
> being taken (without regard for units) from the X11 device.  So
> approximately 7 inches square screen window gets drawn into
> approximately 7 *pixel* square PNG file and (understandably) R complains
> that there is not enough room for the plot.

Yes, that it how it is documented to work.

> Another workaround is something like ...
>
> dev.print(png, width=480, height=480)

(Just one will do if you want to preserve the aspect ratio.)

> ... and a fix requires making dev.print() smarter so that it figures out
> that it needs to convert width/height from inches to pixels.

I don't think there is a way to do that unambiguously (there is no 
standard way to do the conversion), and in any case dev.print() was passed 
a function, not the name of a function, and so does not in general know 
how it behaves (and your 'png' need not be R's png()).

All we can do is to re-emphasize this on the help page, and add a warning 
if a known bitmap device is detected (possibly inaccurately) by name.

BTW, I think it was perverse to call dev.print() except to do printing.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel