Quipu, a simple Nelder Mead solver in F#

Some time back, I wrote a small post digging into the mechanics behind the Nelder Mead solver. As it turns out, I had a use for it recently, and after copy-pasting my own code a few times, I figured it would make my life easier to turn that into a NuGet package, Quipu.

So what does it do, and why might you care?

A code example might be the quickest explanation here. Suppose that, for whatever reason, you were interested in the function f(x) = x ^ 2, and wanted to know for what value of x this function reaches its minimum.

That is easy to solve with the Quipu Nelder-Mead solver:

open Quipu
open Quipu.NelderMead

let f x = x ** 2.0
let solution =
    NelderMead.solve 
        Configuration.defaultValue
        (Objective.from f) [ 100.0 ]
printfn $"{solution}"

… which produces the following result:

Optimal (0.0001556843433, [|0.01247735322|])

The function f reaches a minimum of 0.0001, for x = 0.0124.

Study notes: function minimization with DiffSharp

This post is intended primarily as a note to myself, keeping track as my findings as I dig into automatic differentiation with DiffSharp. Warning: as a result, I won’t make a particular effort at pedagogy – hopefully you’ll still find something of interest in here!

The main question I am interested in here is, how can I use DiffSharp to find the minimum of a function? I will take a look first at basic gradient descent, to get us warmed up. In a future installment I plan to explore using the built-in SGD and Adam optimizers for that same task.

The full code is here on GitHub, available as a .NET interactive notebook.

Test function

The function we will be using in our exploration is the following:

$f(x,y)=0.26 \times (x^2+y^2) - 0.48 \times (x \times y)$

This function, which I lifted this function from this blog post, translates into this F# code:

let f (x: float, y: float) =
    0.26 * (pown x 2 + pown y 2) - 0.48 * x * y

Graphically, this is how the function looks like:

This function has a global minimum for (x = 0.0, y = 0.0), and is unimodal, that is, it has a single peak (or valley in this case). This makes it a good test candidate for function minimization using gradient descent.

Simulating the Wrapinator 5000

It is that time of the year again! The holidays are approaching, and the F# Advent calendar is in full swing. My contribution this year might not be for the broadest audience, sorry about that :) But if you are into F#, probability theory, and numeric optimization, this post is for you - hope you enjoy it! And big shout out to Sergey Tihon for making this happen once again.

You can find the full code for this post here on GitHub.

With the Holidays approaching, Santa Claus, CEO of the Santa Corp, was worried. In preparation for the season’s spike in activity, Santa had invested in the top-of-the-line gift wrapping machine for the Elves factory, the Wrapinator 5000. This beast of a machine has two separate feeders, one delivering Paper, the other Ribbon, allowing the elves to wrap gifts at a cadence never achieved before.

So why worry? Mister Claus, being no fool, had also invested in monitoring, and the logs for the Wrapinator 5000 showed quite a few failures. Would the gift wrapping production lines hold up during the Merry Season?

Mister Claus fiddled anxiously with his luscious beard, alone in his office. And then, being a Man of Science, he did what any self-respecting CEO would do, and decided it was time to build a simulation model of his Wrapinator 5000. With a simulation model in hand, he could analyze his Elves factory, evaluate potential alternative operating policies, and most importantly, get that peace of mind he so desperately longed for.

Maximum Likelihood Estimation of Weibull reliability with DiffSharp

This post is a continuation of my exploration of DiffSharp, an F# Automatic Differentiation library. In the previous post, I covered some introductory elements of Maximum Likelihood Estimation, through a toy problem, estimating the likelihood that a sequence of heads and tails had been generated by a fair coin.

In this post, I will begin diving into the real-world problem that motivated my interest.

The general question I am interested in is about reliability. Imagine that you have a piece of equipment, and that from time to time, a component of that piece of equipment experiences failures. When that happens, you replace the defective component with a new one, restart the equipment, and proceed until the next failure.

The log for such a process would then look something like this:

Jan 07, 2020, 08:37  Component failure
Jan 17, 2020, 13:42  Component failure
Feb 02, 2020, 06:05  Component failure
Feb 06, 2020, 11:17  Component failure

If you wanted to make improvements to the operations of your equipment, one useful technique would be to simulate that equipment and its failures. For instance, you could evaluate a policy such as “we will pre-emptively replace the Component every week, before it fails”. By simulating random possible sequences of failures, you could then evaluate how effective that policy is, compared to, for instance, waiting for the failures to occur naturally.

However, in order to run such a simulation, you would need to have a model for the distributions of these failures.

The log above contains the history of past failures, so, assuming nothing changed in the way the equipment is operated, it is a sample of the failure distribution we are interested in. We should be able to use it, to estimate what parameters are the most likely to have generated that log. Let’s get to it!

Environment: dotnet 6.0.302 / Windows 10, DiffSharp 1.0.7

Full code on GitHub

Weibull time to failure

What we are interested in here is the time to failure, the time our equipment runs from the previous component replacement to its next failure. A classic distribution used to model time to failure is the Weibull distribution.

First look at the new DiffSharp

This post is intended as an exploration of DiffSharp, an Automatic Differentiation, or autodiff F# library. In a nutshell, autodiff allows you to take a function expressed in code - in our case, in F# - and convert it in an F# function that can be differentiated with respect to some parameters.

For a certain niche population, people who care about computing gradients, this is very powerful. Basically, you get gradient descent, a cornerstone of machine learning, for free.

DiffSharp has been around for a while, but has undergone a major overhaul in the recent months. I hadn’t had time to check it out until now, but a project I am currently working on gave me a great excuse to look into these changes. This post will likely expand into more, and is intended as an exploration, trying to figure out how to use the library. As such, my code samples will likely be flawed, and should not be taken as guidance. This is me, getting my hands on a powerful and complex tool and playing with it!

With that disclaimer out of the way, let’s dig into it. In this installment, I will take a toy problem, a much simpler version of the real problem I am after, and try to work my way through it. Using DiffSharp for that specific problem will be total overkill, but will help us introduce some ideas which we will re-use later, for the actual problem I am interested in, estimating the parameters of competing random processes using Maximum Likelihood Estimation.

Environment: dotnet 6.0.302 / Windows 10, DiffSharp 1.0.7

An unfair coin

Let’s consider a coin, which, when tossed in the air, will land sometimes Heads, sometimes Tails. Imagine that our coin is unbalanced, and tends to land more often than not on Heads, 65% of the time:

type Coin = {
    ProbabilityHeads: float
    }

let coin = { ProbabilityHeads = 0.65 }

Supposed now that we tossed that coin in the air 20 times. What should we expect?