# [R] Rmpfr correlation

tri@g m@iii@g oii gvd@et@dk tri@g m@iii@g oii gvd@et@dk
Sun Jul 12 14:57:05 CEST 2020

```Thanks a lot - solved the issue!

BW
Troels

-----Oprindelig meddelelse-----
Sendt: 12. juli 2020 12:59
Til: tring using gvdnet.dk; r-help mailing list <r-help using r-project.org>
Emne: Re: [R] Rmpfr correlation

Hello,

Why not write a function COR? Not one as general purpose as stats::cor but a simple one, to compute the sample Pearson correlation only.

library(Rmpfr)

COR <- function(x, y){
precBits <- getPrec(x)
n <- mpfr(length(x), precBits = precBits)
x.bar <- mean(x)
y.bar <- mean(y)
numer <- sum(x*y) - n*x.bar*y.bar
denom <- sqrt(sum(x*x) - n*x.bar*x.bar) * sqrt(sum(y*y) - n*y.bar*y.bar)
numer/denom
}

set.seed(2020)
KA <- mpfr(10^-4.6, 128)
x <- rnorm(100)*KA
y <- rnorm(100)*x

cor(as.numeric(x), as.numeric(y)) # -0.1874986
# -0.1874986

COR(x, y)
#1 'mpfr' number of precision  128   bits
# -0.1874985950531874160800643775644747505073

Hope this helps,

Às 10:42 de 12/07/20, tring using gvdnet.dk escreveu:
> Dear friends - I'm calculating buffer capacities by different methods and
> need very high precision and package Rmpfr is working beautifully. However,
> I have not been able to find out how to keep precision when finding
> correlations.
>
> library(Rmpfr)
>
> KA <- mpfr(10^-4.6, 128)
>
> x <- rnorm(100)*KA
>
> y <- rnorm(100)*x
>
> cor(x,y) # "x" must be numeric
>
> cor(as.numeric(x),as.numeric(y))# 0.2918954
>
>
>
> In my concrete application I get cor = 1 for
> cor(as.numeric(dff\$BB),as.numeric(BBVS)) even though I have
>
>
>
> str(summary((dff\$BB)-(BBVS)))
> Class 'summaryMpfr' [package "Rmpfr"] of length 6 and precision 128
>   4.61351010833e-8 7.33418976521e-7 1.31009046563e-5 3.76407022709e-5
> 5.72386764888e-5 ...
>
>
>
> I am on windows 10
>
> R version 3.6.1
>
> Best wishes
> Troels Ring,
> Aalborg, Denmark
>
>
>
>
>
>
> This email has been scanned by BullGuard antivirus protection.
> p&url=/>
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help