Darren Wilkinson's blog: Recent Episodes

None

Statistics, computing, functional programming, data science, Bayes, stochastic modelling, systems biology and bioinformatics

View Details

Breeze is the standard scientific and numerical library for Scala. For linear algebra operations, it builds on top of the Java library, netlib. This provides a nice interface to BLAS and related libraries which allows the use of native optimised libraries and will also gracefully fall back to using pure Java implementations if optimised native code libraries can’t be found. This is great, since it leads to good code portability, but the Java implementations will typically be slower than optimised native libraries for large matrices, so if you care about speed it is important to install optimised libraries on your system and configure netlib to use them.

See the netlib Readme for details of installing native libraries and setting the relevant system properties. In many cases netlib will automatically detect and use native libraries, but it’s not foolproof, especially if you are using nonstandard libraries, or if you have multiple libraries installed and you want to specify which one to use. Briefly, you can override default settings by setting the properties blas, lapack and arpack. Each of this can be set by using either nativeLib, to specify the name of a library in your system library search path or nativeLibPath, to set the full path to the library you require. Full examples of the two approaches are:

-Ddev.ludovic.netlib.blas.nativeLib=libopenblas.so-Ddev.ludovic.netlib.blas.nativeLibPath=/usr/lib/x86_64-linux-gnu/libopenblas.so Obviously these need to be customised to your requirements. lapack and arpack properties are set similarly.

What the netlib readme doesn’t discuss is how to set these properties in Scala projects, or how to check/verify the libraries being used in Scala Breeze projects. These are discussed below.

Setting netlib properties in Scala projectsExactly how you set system properties depends on exactly how you are building and running your Scala code. However, some build tools will read the environment variable JAVA_OPTS, so setting this will very often work. For example, if you are using a bash-like shell, you could set a relevant property with something like:

export JAVA_OPTS="-Ddev.ludovic.netlib.blas.nativeLib=libblas.so" Multiple properties should be separated with a space.

export JAVA_OPTS="-Ddev.ludovic.netlib.blas.nativeLib=libblas.so -Ddev.ludovic.netlib.lapack.nativeLib=liblapack.so" Then this might be picked up and used by your build tool. If so, this is often the preferred approach.

scala-cliscala-cli is a popular tool for compiling and running small Scala projects. You can pass in a property directly at the end of the scala-cli command-line:

scala-cli run breeze-test.scala '-Ddev.ludovic.netlib.blas.nativeLib=libblas.so' Multiple properties should be included separately:

scala-cli run breeze-test.scala '-Ddev.ludovic.netlib.blas.nativeLib=libblas.so' '-Ddev.ludovic.netlib.lapack.nativeLib=liblapack.so' If you prefer, you can include the option in your Scala source code in the scala-cli headers, which might then look similar to:

//> using scala 3.3.0//> using dep org.scalanlp::breeze:2.1.0//> using dep org.scalanlp::breeze-viz:2.1.0//> using javaOpt -Ddev.ludovic.netlib.blas.nativeLibPath=/usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0 The (significant) disadvantage of this approach is that it makes the code less portable.

sbtMany larger Scala project are built using sbt.sbt checks the JAVA_OPTS environment variable, so this is often the preferred way to configure the libraries that you want to use.

Alternatively, you can explicitly set them at the sbt command line by inserting them before the required task:

sbt -Ddev.ludovic.netlib.blas.nativeLib=libblas.so run If you prefer, you can include the options inside your build.sbt file with something like:

javaOptions ++= Seq( "-Ddev.ludovic.netlib.blas.nativeLibPath=/usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0")fork := true Note that you will probably need to fork the project, as the options will be applied to the forked process. But in addition to the usual pros and cons of forking sbt projects, this approach has the disadvantage of making the code less portable.

Verifying netlib instances in Scala Breeze projectsSince the linking to specific libraries happens at runtime, it is often desirable to be able to check whether native libraries are being used within a running Scala Breeze application. There are various ways to do this, but the basic idea can be illustrated with the following Scala code snippet.

import dev.ludovic.netlib.blas.BLASprint("BLAS: ")println(BLAS.getInstance().getClass().getName())import dev.ludovic.netlib.lapack.LAPACKprint("LAPACK: ")println(LAPACK.getInstance().getClass().getName())import dev.ludovic.netlib.arpack.ARPACKprint("ARPACK: ")println(ARPACK.getInstance().getClass().getName()) If the output from this snippet is something like:

BLAS: dev.ludovic.netlib.blas.JNIBLASLAPACK: dev.ludovic.netlib.lapack.JNILAPACKARPACK: dev.ludovic.netlib.arpack.JNIARPACK then you are using native libraries. The letters JNI in the final word indicate the use of the “Java native interface”. Any other output indicates that a relevant native library has not been found, and the precise output will give some indication of exactly what pure Java implementation has been fallen back to.

If you wanted some more friendly output, you can do something like:

val blas = BLAS.getInstance().getClass().getName() println(s"Using BLAS: $blas") blas match case "dev.ludovic.netlib.blas.JNIBLAS" => println("This is a native BLAS of some sort") case "dev.ludovic.netlib.blas.VectorBLAS" => println("This is the VectorBLAS for Java 16+") case _ => println("Fallen back to a Java BLAS of some sort (probably slow)") Note that this also detects the use of VectorBLAS, discussed below.

Java 16+As explained in the netlib readme, there is an additional wrinkle if you are using a recent JVM (version 16 or higher). Since recent JVMs expose vector operations, it is now possible to write pure Java BLAS libraries with similar performance to native libraries. Recent versions of netlib include such an implementation, called VectorBLAS. So, if you run your Scala Breeze application on a recent JVM, netlib will first check to see if it can detect a VectorBLAS, and if it finds it, it will use it in preference to a native library. However, if it can’t (which is likely to be the case by default), then it will next look for a native BLAS, before eventually falling back to a less performant Java BLAS library. So if you have a good native BLAS installed and configured, then you are probably happy with this, and can safely ignore any errors or warnings about not being able to find a VectorBLAS.

View Details

Functional programming (FP) languages are great for statistical computing, computational statistics, and machine learning. They are particularly well-suited to scalable computation, where this could either mean scaling up to distributed algorithms for very big data, or running algorithms for more moderately sized data sets very fast in parallel on GPUs. However, people unfamiliar with FP often find FP languages quite intimidating, due to the fairly steep initial learning curve. This issue is exacerbated by the fact that there is very little introductory documentation available for people new to FP who are interested in applications to statistical computing and machine learning (ML).

So for some time I’ve been meaning to put together materials for a short course (suitable for self-study) that will get people started with FP in few different languages, with a very basic problem from statistical computing used as the running example, together with a catalogue of resources for further learning, in order to provide people with the information they need to keep going once they have got over the initial hurdle. But as with many things, it never got high enough up my priority list to actually sit down and do it. Fortunately, StatML invited me to deliver some training in advanced statistical computing, so this gave me the perfect motivation to actually assemble something. The in-person training has been delayed (due to the UCU strike), but the materials are all prepared and publicly available, and suitable for self-study now.

The course gives a very quick introduction to the ideas of FP, followed by very quick hands-on introductions to my favourite FP languages/libraries: Scala, Haskell, JAX and Dex. There is also a brief introduction to splittable random number generators which are becoming increasingly popular for the development of functional parallel Monte Carlo algorithms.

If you’ve been vaguely interested in FP for statistical computing and ML but not sure how to get started, hopefully this solves the problem.

An introduction to functional programming for scalable statistical computing and machine learning (short course)

View Details

In June this year the (twice COVID-delayed) Richard J Boys Memorial Workshop finally took place, celebrating the life and work of my former colleague and collaborator, who died suddenly in 2019 (obituary). I completed the programme of talks by delivering the inaugural RSS North East Richard Boys lecture. For this, I decided that it would … Continue reading MCMC code for Bayesian inference for a discretely observed stochastic kinetic model

View Details

Yesterday there was an RSS Read Paper meeting for the paper Unbiased Markov chain Monte Carlo with couplings by Pierre Jacob, John O’Leary and Yves F. Atchadé. The paper addresses the bias in MCMC estimates due to lack of convergence to equilibrium (the “burn-in” problem), and shows how it is possible to modify MCMC algorithms … Continue reading Unbiased MCMC with couplings

View Details

This is the 75th post to this blog. Every 25 posts I produce an index of posts so far for easy reference. If I make it to post 100 I’ll do something similar. 25. Catalogue of my first 25 blog posts 50. Index to first 50 posts 51. Calling Scala code from R using rscala … Continue reading Index to first 75 posts

View Details

Introduction In the previous post I showed how to write your own general-purpose monadic probabilistic programming language from scratch in 50 lines of (Scala) code. That post is a pre-requisite for this one, so if you haven’t read it, go back and have a quick skim through it before proceeding. In that post I tried … Continue reading A probability monad for the bootstrap particle filter

View Details

Background In May I attended a great workshop on advances and challenges in machine learning languages at the CMS in Cambridge. There was an a good mix of people from different disciplines, and a bit of a theme around probabilistic programming. The workshop schedule includes links to many of the presentations, and is generally worth … Continue reading Write your own general-purpose monadic probabilistic programming language from scratch in 50 lines of (Scala) code

View Details

Introduction There is a fairly large literature on reaction-diffusion modelling using partial differential equations (PDEs). There is also a fairly large literature on stochastic modelling of coupled chemical reactions, which account for the discreteness of reacting species at low concentrations. There is some literature on combining the two, to form stochastic reaction-diffusion systems, but much … Continue reading Stochastic reaction-diffusion modelling

View Details

In the previous post I gave a very quick introduction to the smfsb R package. As mentioned in that post, although good for teaching and learning, R isn’t a great language for serious scientific computing or computational statistics. So for the publication of the third edition of my textbook, Stochastic modelling for systems biology, I … Continue reading The scala-smfsb library

View Details

Introduction In the previous post I gave a brief introduction to the third edition of my textbook, Stochastic modelling for systems biology. The algorithms described in the book are illustrated by implementations in R. These implementations are collected together in an R package on CRAN called smfsb. This post will provide a brief introduction to … Continue reading The smfsb R package

View Details

The third edition of my textbook, Stochastic Modelling for Systems Biology has recently been published by Chapman & Hall/CRC Press. The book has ISBN-10 113854928-2 and ISBN-13 978-113854928-9. It can be ordered from CRC Press, Amazon.com, Amazon.co.uk and similar book sellers. I was fairly happy with the way that the second edition, published in 2011, … Continue reading Stochastic Modelling for Systems Biology, third edition

View Details

Introduction In the previous post I gave a brief introduction to Rainier, a new HMC-based probabilistic programming library/DSL for Scala. In that post I assumed that people were using the latest source version of the library. Since then, version 0.1.1 of the library has been released, so in this post I will demonstrate use of … Continue reading Bayesian hierarchical modelling with Rainier

View Details

Introduction Rainier is an interesting new probabilistic programming library for Scala recently open-sourced by Stripe. Probabilistic programming languages provide a computational framework for building and fitting Bayesian models to data. There are many interesting probabilistic programming languages, and there is currently a lot of interesting innovation happening with probabilistic programming languages embedded in strongly typed … Continue reading Monadic probabilistic programming in Scala with Rainier