This is a sample Sage notebook showing a simple linear regression.

In a Sage Notebook, a Python computation is structured as a set of individually executable blocks of code. The state of memory is maintained across the entire running life of the notebook, so code blocks can see the results of previously evaluated blocks in the same notebook.

My first code block just imports the Python Random library for random number generation.

︡ea7a639b-e7c9-4dcb-8d44-b70d0734a278︡{"html":"This is a sample Sage notebook showing a simple linear regression.

\nIn a Sage Notebook, a Python computation is structured as a set of individually executable blocks of code. The state of memory is maintained across the entire running life of the notebook, so code blocks can see the results of previously evaluated blocks in the same notebook.

\nMy first code block just imports the Python Random library for random number generation.

\n\n"}︡ ︠c330de5c-177c-48d2-a836-e348deac4ed6s︠ import random ︡576cc7dc-32e4-4fea-aeb0-0c73a665e421︡ ︠b0b3c49a-6e7d-4644-889f-6b9999341bc4si︠ %htmlThe next block defines a linear relationship with randomly chosen slope and intercept parameters. Note that the slope and intercept change each time this block is evaluated. Once this block is evaluated, however, the function y(x) is fixed for the rest of the notebook.

︡e3720596-98d8-435a-a3ab-c4fa868f0247︡{"html":"The next block defines a linear relationship with randomly chosen slope and intercept parameters. Note that the slope and intercept change each time this block is evaluated. Once this block is evaluated, however, the function y(x) is fixed for the rest of the notebook.

\n\n"}︡ ︠4bdb521d-1344-48a7-b944-22b54ee231b9s︠ intercept = random.uniform(-50, 50) slope = random.uniform(-2, 2) def y(x): return slope*x + intercept print [slope, intercept] ︡392a05a4-8dd5-4911-b9a7-e524c26aa8aa︡{"stdout":"[-0.7252207008692442, 16.54838088794004]\n"}︡ ︠dfb5519b-1e6c-45ea-99a4-818857945ecbsi︠ %htmlCreate some noisy data where the X values are chosen from a normal distribution and the Y values are y(x) plus a normally-distributed noise term.

︡bf972800-92b4-4a14-9031-693fe8ed2259︡{"html":"Create some noisy data where the X values are chosen from a normal distribution and the Y values are y(x) plus a normally-distributed noise term.

\n\n"}︡ ︠ea9c1210-0b1d-4ac2-a33d-94e49f6ec50cs︠ X = [random.normalvariate(100, 20) for i in xrange(1000)] Y = [y(x) + random.normalvariate(0, 10) for x in X] ︡0c2fb9db-523b-42d7-bbcd-0373a524fd33︡ ︠34a93121-1a05-4629-89d8-ba63093ea5dbsi︠ %htmlIt's easy to plot the data right here in the notebook. Note for future reference, the code

zip(X, Y)

is a built-in Python function that takes two lists X=[x1, x2, ...] and Y=[y1, y2, ...] and munges them into a single list of tuples [(x1, y1), (x2, y2), ...]

︡28f78a1b-5e6e-4ac1-9441-a2787cb7f690︡{"html":"It's easy to plot the data right here in the notebook. Note for future reference, the code

\nzip(X, Y)\n

is a built-in Python function that takes two lists X=[x1, x2, ...] and Y=[y1, y2, ...] and munges them into a single list of tuples [(x1, y1), (x2, y2), ...]

\n\n"}︡ ︠4a1a9089-4c4c-4cd3-8c88-c68a8d75c007s︠ list_plot(zip(X,Y)) ︡945516a7-4801-4767-9b77-74a3ac8ac43b︡{"once":false,"file":{"show":true,"uuid":"d77553bd-fd4a-457c-8663-b63aa2ccb584","filename":"/projects/14f921e7-850e-496d-a966-c5f281a9e725/.sage/temp/compute2-us/26744/tmp_kryO3C.svg"}}︡{"html":""}︡ ︠6e6f357e-8061-4957-bf58-65b005649f01si︠ %htmlSage comes pre-loaded with tons of Python libraries for scientific computing and applied math. Here, I'm just using a simple function called find_fit to fit a model (linear model in this case) to my data. Compare the a and b parameters of my model with the slope and intercept printed above.

︡b6c9adf3-8d08-4b83-af81-b40253b607a1︡{"html":"Sage comes pre-loaded with tons of Python libraries for scientific computing and applied math. Here, I'm just using a simple function called find_fit to fit a model (linear model in this case) to my data. Compare the a and b parameters of my model with the slope and intercept printed above.

\n\n"}︡ ︠de179cee-294e-474b-84ab-d7df31d4001ds︠ var('a, b') model(i) = a*i + b find_fit(zip(X,Y), model) ︡50ef6964-edfa-4356-88df-1735fa72358b︡{"stdout":"(a, b)\n"}︡{"stdout":"[a == -0.7121013864723931, b == 15.01876764562296]\n"}︡ ︠08091796-6d11-414b-ba2c-a8efb59b56c3si︠ %htmlShow the fit model and my data together in a plot.

︡e84bab46-12fa-468e-b75b-e17082020ae7︡{"html":"Show the fit model and my data together in a plot.

\n\n"}︡ ︠053a4e20-10d7-4cda-861a-5f883eb0c1fas︠ points(zip(X,Y))+plot(model(a=find_fit(zip(X,Y),model)[0].rhs(),b=find_fit(zip(X,Y),model)[1].rhs()), xmin=min(X), xmax=max(X),color='red') ︡b770514f-6b89-475e-9fb7-71c4a95846cf︡{"once":false,"file":{"show":true,"uuid":"dc220646-4997-4af2-91d0-1a46a0647d2f","filename":"/projects/14f921e7-850e-496d-a966-c5f281a9e725/.sage/temp/compute2-us/26744/tmp_Ppzoto.svg"}}︡{"html":""}︡ ︠734727bd-1c17-4a08-9dfd-c6b1b970cc69︠