Evolving notes, images and sounds by Luis Apiolaza

Category: r (Page 2 of 20)

Creating an n x n autocorrelation matrix

Between covid-19 news and announcements of imminent Russia-Ukraine wars I needed a bit of a distraction. Sooo, here it is how to create an n x n autocorrelation matrix based on a correlation rho, with a simple 5 x 5 example in R:

n <- 5
rho <- 0.9
mat <- diag(n)
mat <- rho^abs(row(mat) - col(mat))

This produces the following output:

       [,1]  [,2] [,3]  [,4]   [,5]
[1,] 1.0000 0.900 0.81 0.729 0.6561
[2,] 0.9000 1.000 0.90 0.810 0.7290
[3,] 0.8100 0.900 1.00 0.900 0.8100
[4,] 0.7290 0.810 0.90 1.000 0.9000
[5,] 0.6561 0.729 0.81 0.900 1.0000

How does this work? Starting from defining n as 5, and the basic rho correlation as 0.9, then mat <- diag(n) creates a diagonal matrix—a square matrix with ones in the diagonal and zeros elsewhere—with 5 rows and 5 columns.

Continue reading

The beauty of code vectorisation

I came across this problem in Twitter:

Description of the problem by Fermat’s Library.

The basic pieces of the problem are:

  • We need to generate pseudorandom numbers with an identical distribution, add them up until they go over 1 and report back how many numbers we needed to achieve that condition.
  • Well, do the above “quite a few times” and take the average, which should converge to the number e (2.7182…).
Continue reading

Recreational programming

I think programming, aka coding, is a fun activity. We are solving a problem subject to a set of constraints that can be time, memory, quantity of code, language, etc. Besides being a part of my work, coding is also a good distraction when doing it for the sake of it.

In this type of situation I like to set myself additional constraints. For example, I try to use only base R or a very minimal set of features. Moreover, I prefer vectorised code: there are no explicit loops in the code. For toy problems this restriction doesn’t make much of a difference in execution time, but in real problems R (and Matlab) run much more slowly with non-vectorised code.

Continue reading

Fixing Rcpp warning in Mac OS

In Mac OS I was getting an annoying warning when compiling Cpp code via Rcpp in R:

ld: warning: directory not found for option '-L/usr/local/gfortran/lib/gcc/x86_64-apple-darwin15/version_number'

Adding a file called Makevars in ~/.R defining FLIBS to the actual location of gfortran in my machine fixed the problem. In my case, Makevars only contains the following line:

FLIBS=-L/usr/local/gfortran/lib

That’s it.

Implementing a model as an R package

In our research group we often have people creating statistical models that end up in publications but, most of the time, the practical implementation of those models is lacking. I mean, we have a bunch of barely functioning code that is very difficult to use in a reliable way in operations of the breeding programs. I was very keen on continue using one of the models in our research, enough to rewrite and document the model fitting, and then create another package for using the.model in operations.

Unfortunately, neither the data nor the model are mine to give away, so I can’t share them (yet). But I hope these notes will help you in you are in the same boat and need to use your models (or ‘you’ are in fact future me, who tend to forget how or why I wrote code in a specific way).

Continue reading
« Older posts Newer posts »

© 2024 Palimpsest

Theme by Anders NorenUp ↑