In [1]: 
using PyPlot

Loading help data...
In [31]: 
# a single trial of 100 data points
X = rand(100)
E = randn(100) * 0.5
Y = 2 + 3X + E

A = hcat(ones(100), X)
U, W, V = svd(A)
a = V * diagm(1 ./ W) * U' * Y

scatter(X, Y)
plot_x = linspace(0, 1, 100)
plot(plot_x, a[1] + a[2] * plot_x)
title("a[1] = $(a[1])\na[2] = $(a[2])")
Out[31]:
PyObject <matplotlib.text.Text object at 0x11df31b90>
In [78]: 
# 100 bootstrapped data sets
a = zeros(2)
plot_x = linspace(0, 1, 100)

X = rand(100)
E = randn(100) * 0.5
Y = 2 + 3X + E

for i in 1:100

    idxs = rand(1:100, 100)
    X_ = X[idxs]
    Y_ = Y[idxs]
    A = hcat(ones(100), X_)
    U, W, V = svd(A)
    a_ = V * diagm(1 ./ W) * U' * Y_
    plot(plot_x, a_[1] + a_[2] * plot_x, c="black", alpha=0.3)
    scatter(X_, Y_, alpha=0.05)
    a += a_
end
a /= 100


plot(plot_x, a[1] + a[2] * plot_x)
title("a[1] = $(a[1])\na[2] = $(a[2])")
Out[78]:
PyObject <matplotlib.text.Text object at 0x1246d9590>
In [48]: 
# 100 independent data sets
a = zeros(2)
plot_x = linspace(0, 1, 100)
for i in 1:100
    X = rand(100)
    E = randn(100) * 0.5
    Y = 2 + 3X + E

    A = hcat(ones(100), X)
    U, W, V = svd(A)
    a_ = V * diagm(1 ./ W) * U' * Y
    plot(plot_x, a_[1] + a_[2] * plot_x, c="black", alpha=0.3)
    scatter(X, Y, alpha=0.05)
    a += a_
end
a /= 100


plot(plot_x, a[1] + a[2] * plot_x)
title("a[1] = $(a[1])\na[2] = $(a[2])")
Out[48]:
PyObject <matplotlib.text.Text object at 0x122db5b90>
In [27]: 
# Trying out a quadratic fit for fun
X = rand(100)
E = randn(100) * 0.5
Y = 2 + 3X + E

A = hcat(ones(100), X, X .* X)
U, W, V = svd(A)
a = V * diagm(1 ./ W) * U' * Y

scatter(X, Y)
plot_x = linspace(0, 1, 100)
plot(plot_x, a[1] + a[2] * plot_x + a[3] * plot_x .* plot_x)
title("a[1] = $(a[1])\na[2] = $(a[2])\na[3] = $(a[3])")
Out[27]:
PyObject <matplotlib.text.Text object at 0x11d7d2110>