gistfile1.txt

#r @"..\packages\R.NET.1.5.5\lib\net40\RDotNet.dll"
#r @"..\packages\RDotNet.FSharp.0.1.2.1\lib\net40\RDotNet.FSharp.dll"
#r @"..\packages\R.NET.1.5.5\lib\net40\RDotNet.NativeLibrary.dll"
#r @"..\packages\RProvider.1.0.2\lib\RProvider.dll"

open System
open RDotNet
open RProvider

open RProvider.``base``
open RProvider.graphics
open RProvider.stats

let rng = Random()
let rand () = rng.NextDouble()

// Create artificial data:
// the model formulation is
// Y = 5.0 + 3.0 * X1 - 2.0 * X2 + Epsilon (noise)
let X1s = [ for i in 0 .. 9 -> 10. * rand () ]
let X2s = [ for i in 0 .. 9 -> 5. * rand () ]
let Ys = [ for i in 0 .. 9 -> 5. + 3. * X1s.[i] - 2. * X2s.[i] + rand () ]

// Create a dataset
let dataset =
    namedParams [
        "Y", box Ys;
        "X1", box X1s;
        "X2", box X2s; ]
    |> R.data_frame

// Run a linear regression
// See R docs for the somewhat esoteric construction of formulas:
// http://stat.ethz.ch/R-manual/R-patched/library/stats/html/formula.html
let result = R.lm(formula = "Y~X1+X2", data = dataset)

// Retrieve the coefficients and residuals:
// The coefficients should be around 5.0, 3.0 and -2.0
// According to R docs, result is a list;
// We can retrieve elements by their name
// http://stat.ethz.ch/R-manual/R-patched/library/stats/html/lm.html
let coefficients = result.AsList().["coefficients"].AsNumeric()
let residuals = result.AsList().["residuals"].AsNumeric()

// We can also get some summary statistics on our model,
// like R2, which measures goodness-of-fit.
// See R docs for the details on Summary:
// http://stat.ethz.ch/R-manual/R-patched/library/stats/html/summary.lm.html
let summary = R.summary(result)
summary.AsList().["r.squared"].AsNumeric()

// We can also pass results to R.plot - fancy!
R.plot result