[R] Problem with lm Giving Wrong Results

[Labone, Thomas]

Thanks.

Can someone point me to the procedure for switching from the Intel Math Library back to the standard math library so that I can see if the problem is associated with using MKL?


Thomas R. LaBone
PhD student
Department of Epidemiology and Biostatistics
Arnold School of Public Health
University of South Carolina
Columbia, South Carolina USA


________________________________
From: J C Nash <profjcnash using gmail.com>
Sent: Thursday, December 2, 2021 9:31 AM
To: Labone, Thomas <labone using email.sc.edu>; r-help using r-project.org <r-help using r-project.org>
Subject: Re: [R] Problem with lm Giving Wrong Results

I get two similar graphs.

https://protect2.fireeye.com/v1/url?k=99a20372-c6393a53-99a24db3-862c53b6784d-c549bca0bcc4210d&q=1&e=8ce563a7-fac1-41d3-8699-13c4a9417e00&u=https%3A%2F%2Fweb.ncf.ca%2Fnashjc%2Fjfiles%2FRplot-Labone-4095.pdf
https://protect2.fireeye.com/v1/url?k=6d56bd36-32cd8417-6d56f3f7-862c53b6784d-85baef704a412b23&q=1&e=8ce563a7-fac1-41d3-8699-13c4a9417e00&u=https%3A%2F%2Fweb.ncf.ca%2Fnashjc%2Fjfiles%2FRplotLabone10K.pdf

Context:

R version 4.1.2 (2021-11-01)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 20.2

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/liblapack.so.3

locale:
  [1] LC_CTYPE=en_CA.UTF-8       LC_NUMERIC=C               LC_TIME=en_CA.UTF-8        LC_COLLATE=en_CA.UTF-8
  [5] LC_MONETARY=en_CA.UTF-8    LC_MESSAGES=en_CA.UTF-8    LC_PAPER=en_CA.UTF-8       LC_NAME=C
  [9] LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C

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

loaded via a namespace (and not attached):
  [1] compiler_4.1.2  fastmap_1.1.0   htmltools_0.5.2 tools_4.1.2     yaml_2.2.1      rmarkdown_2.11  knitr_1.36
  [8] xfun_0.28       digest_0.6.28   rlang_0.4.12    evaluate_0.14
 >

Hope this helps,

JN


On 2021-12-02 5:50 a.m., Labone, Thomas wrote:
> #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> # This works
> n <- 1000# OK <= 4095
> Z <- qnorm(ppoints(n))
>
> k <- sort(rlnorm(n,log(2131),log(1.61)) / rlnorm(n,log(355),log(1.61)))
>
> quantile(k,probs=c(0.025,0.5,0.975))
> summary(k)
>
> fit <- lm(log(k) ~ Z)
> summary(fit)
>
> gm <- exp(coef(fit)[1])
> gsd <- exp(coef(fit)[2])
> gm
> gsd
>
> plot(Z,k,log="y",xlim=c(-4,4),ylim=c(0.1,100))
> lines(Z,gm*gsd^Z,col="red")
>
> #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> #this does not
> n <- 10000# fails >= 4096 = 2^12
> Z <- qnorm(ppoints(n))
>
> k <- sort(rlnorm(n,log(2131),log(1.61)) / rlnorm(n,log(355),log(1.61)))
>
> quantile(k,probs=c(0.025,0.5,0.975))
> summary(k)
>
> fit <- lm(log(k) ~ Z)
> summary(fit)
>
> gm <- exp(coef(fit)[1])
> gsd <- exp(coef(fit)[2])
> gm
> gsd
>
> plot(Z,k,log="y",xlim=c(-4,4),ylim=c(0.1,100))
> lines(Z,gm*gsd^Z,col="red")
>
>
> #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

	[[alternative HTML version deleted]]