Skip to content

Instantly share code, notes, and snippets.

@benmoran
Created April 28, 2014 21:10
Show Gist options
  • Star 1 You must be signed in to star a gist
  • Fork 1 You must be signed in to fork a gist
  • Save benmoran/49eadfdfab4fdfc3ef53 to your computer and use it in GitHub Desktop.
Save benmoran/49eadfdfab4fdfc3ef53 to your computer and use it in GitHub Desktop.
Sympy code gen demo
{
"metadata": {
"name": ""
},
"nbformat": 3,
"nbformat_minor": 0,
"worksheets": [
{
"cells": [
{
"cell_type": "code",
"collapsed": false,
"input": [
"import sympy as S\n",
"S.var(\"x mu sigma\", real=True)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\begin{pmatrix}x, & \\mu, & \\sigma\\end{pmatrix}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAGgAAAAaBAMAAACp2ARYAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMkS7zRCZdiKJ71Rm\nq90icBAQAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABnUlEQVQ4EZVSMUvDQBT+0jZt2oYojh3UTRex\n4OIoCG5CBhdBaBGl6GKXWioIFX9Au4ibdnNwySIIgmYWFUFBwaWKutoWq0KG+i5pmkum3C33ve99\n37177w7ACIRWkamlrJAHU0we08VMSZP0BTEP8AEoXVHTNpDMi5pqQMoQNWkmxgTnAMSzuKZCcvl5\n4jBQMAPsBih5s9f7ARINrFFmGXfF94DiCPgMUKOrl7NVINLFDGUOcKyfBxQLkC0/JZ9CbREltbFI\nm455f55ubCEdeIz4HCtCj9Tp3+E3aCIHiXxryESCVVIsuxKrKPsEiDeh1as+rjaM2D0zdWzTGZ1Q\n8gmgZVEz6j4upyOjEyP94YHGYUVbioncHKfJNbG1Z8jfHBU1lDcWUmM71PRTqXIDaHzjVxvlBD3d\nFzu6v+RCxY7UBqY9etJN0/7o4EiV4/qQfkRqfECbAwScOFjlKBfGiojX3UA2XUQTajt4xaMGiAYi\nNdyIPzXSdFjTTXL7LeEXLg4FX0m1FErpiVTWTzrvEWHQvi26CCP1NOseFEL/hdpXxpw09kcAAAAA\nSUVORK5CYII=\n",
"prompt_number": 1,
"text": [
"(x, \u03bc, \u03c3)"
]
}
],
"prompt_number": 1
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here's the likelihood for a 1-d normal distribution, in terms of two parameters: the mean $\\mu$ and standard deviation $\\sigma$"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"\n",
"f = 1/(sqrt(2*S.pi)*sigma)*S.exp(-(x-mu)**2/(2*sigma**2))\n",
"f\n"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\frac{\\sqrt{2}}{2 \\sqrt{\\pi} \\sigma} e^{- \\frac{\\left(- \\mu + x\\right)^{2}}{2 \\sigma^{2}}}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAIsAAAA0BAMAAABIq9ZsAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAInarRM2ZVBDdiWbv\nuzJCz3LGAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADJUlEQVRIDaVWS2gTQRj+kuwmm7apsYKICH2I\nIK1IoILPQy4eRMQIHjx4SKu0B197UehBXJQK9mIED1IQgzcrYqg3vXgTPLXgVSkIFUWhoEgtaPxn\ndyY7s7ObZjdzyH7f9z92ZnYy/w8oY6DJxpqiJSAHEsToIamKriVQsnaCID3kni5tqvQvLAd9poHt\nVz4F1fa8gBMBB6sGq4SX5YDcnpq476gePQ1kHPS9VdVN2UzAg7amr4bMv4Ds0dyKLPdwMgxYZ2QD\n4THAXI9Kc1JxznKWAe44vsUoE665PL/uyxIqSRgQaYyGNTHnW/bRwUsPufx03Zd9lKpynC4XCLE0\nA09fAUv55h9uAVK3fwF5jz5uqTIwhjgziz+QHT88fhD2Kdoher00rL/Aa5dnlyTZh702rKmpqWn0\nYZFUd1EjBNQ0uFbEIzdotx8qI6PB2XGctb009Q85WpTsBAxWe+tMyVWxVbV4jM6lN56NjJa9NEdH\n9wNC5tbMet5hcBWgsxAyhL93SsSXypVV1/6NG0xIzS9cqKkWzgb4s+I+Lc5EOk6B8+cYNOnyC09j\nlFuuEtgrYRfeqgYVj9/8SMOBe7tG/az5oWbRx52gJSN81p3ESj4VK/xPI7m0g99e8GnnuppNazGz\nHe/C3MiENjGxGOuQZooQjAfYqZtoMdeb87i7o67bQpVdtnVVN8wWcyW6DJ4o5yvqKLASfGxxzNHS\n0GJMx9R1zbEl/GwhCdBiTGyThE3hRpgHLcZ4Xg+zRGnv6WKNsrXTB72dEi57YOkfXBjbPAN9SXry\nTRtn1SQV7S76Erlod9GXyEWbim/SIRft6aRJlKJNfUk3QxRt6kvS7ItXYV1uNn/HTSmKNm3Nwpfh\nRSpTXy/OLHd8m4jXiaI9hn67p0glzHqH3jVh7fTJirZbPdjWFEC3UbYCI/ZFy4p2qy/5zJqcLbRJ\ncWfDirbfl1wCTWNwBWap08Vwv1VWtEVfQkWDGh1qLb878dLwos37knQNdDsU7NTDeFlE0eZ9SaGB\nI/SlJkdjTka8VPQlgid88r4kYbQf5vUlPk+IovqSmOni9iVq+v8RVuYGHnhB5AAAAABJRU5ErkJg\ngg==\n",
"prompt_number": 2,
"text": [
" 2 \n",
" -(-\u03bc + x) \n",
" \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 2 \n",
" ___ 2\u22c5\u03c3 \n",
"\u2572\u2571 2 \u22c5\u212f \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" ___ \n",
" 2\u22c5\u2572\u2571 \u03c0 \u22c5\u03c3 "
]
}
],
"prompt_number": 2
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And here's the log likelihood, which we can tidy up:"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll = log(f)\n",
"ll"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\log{\\left (\\frac{\\sqrt{2}}{2 \\sqrt{\\pi} \\sigma} e^{- \\frac{\\left(- \\mu + x\\right)^{2}}{2 \\sigma^{2}}} \\right )}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAMsAAAA/BAMAAAC81iEPAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMqvNiRDduyJ2RFSZ\n72bxr6VbAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAFqUlEQVRYCa1XXWgcVRQ+s7OzOzu7ycbgmw9Z\nihXpS4J/BKXsPkSkGpo1loCQkqVYtCAmCBLFBxelpWAlEQIaqcngQ0EL7YpQQQoZ7EOFBhtpLYVi\nsw/+EHxIUhPbxtT13Hvn587cO5OZtRd25pzv+849s3fu3HsugLStSNE2QMOKCMrXIshk1DsR8o+6\nIshk1EEzXH8knErKGNdDI9J9AjXVIk2AYwAXQjVTpkDtEZC4wFgjTDkpELr4/wRNCJBeCiHy6wJR\nCH0kQeoCysDbxNbuuojfUMVHHwF9cbLLL9vJy8IpKhmqypU9TQG/BMOQ3hLgSECFMZpgSuyOxvUK\n0VoNdgP8IuA7AOcoLxkcgmu3hOh0E64B9MpH7WufvMK8XB076qd27g6DAtdOcXBGAMrVkDRG0xdv\n2V4JYJi9FG3DJ3Cc7Jpjufc5ao2zMBdlRsrvW7Z7BrQfDjB7r19he8U+HjZM9GoEMeRPdYJwpH0G\nE3i1APSZuRJkIdP6m8AAQ13s7r+erfD+55g0N02QDOlFbF86UL8yAd/vOrSrkTcmKqA6MMACDfd8\nZg2ZHKK/jG8qQ4GLHMyZfwKMLi4uNmCpk3Rn4Y88FpdmTPp8+LK5pt0G+JX4uRKHcuYZ2zb61NNo\nWgDprNnAQXNbx5JrcoY9og5yuQ6HiX0e9LqD8fdh20kfmX0eTQsgtbKvC7iZkZLNaB0fn29TVmcF\nfcUC8i2IzelPZaxlK77zlOl/Pdu18tuuSY3UZoaM4rGBwUf8hO3p08ywX4Yz4pYn7pRN0VxgRVU2\nzpOIcsuZn148s74KAsTPOenQNnBZ6X4ioCoE/+LeJwMKn1uk22rYhUr1v3DhOeSLAsgG15oFK6Bg\nbn4S21siZW8zPHEPnTd5AG01OC/sdxuQhbrONsMJSA0RTNOxyQnaMJ1thgtdxhcVTFNc4wRJzAef\nqzP5uWDUOBKY5gD5sL4YXPmE8HHSaPtOYmSgTRg1itjbDMeyNMprkJ02NuGNHKHipDlV1ytcN8zs\n09hw29sMx5e7yL/JVEBZT9WAfX5jOw+afgNGMTLYsPL+pvWKu814dNkkaXoaABuFJSfNEuUjvgdI\n3Zz9zevEtV6q50sZb5txcXjPJGkW8HpbX9duUGKMpfFEolUsiRh+gbtBrapVCcXSkH+zrc/MdFFF\njEErViR9wejRigrvyxiWpsMCYytv2oIYU6CjAjAt9Dfeqhn7kRFb2QSYBOUnUCvG0w/UqSBGmlwf\njDIxjdjxgjOtu/UonBjYQ5bgexYJkKXR59na6HR4cq7pmHHu4+x1EGm2oX1IT0/OYsMXzSkzTm+h\nGvJ52u1bvP9B7Iy9pvFF8we2qM3bshc3gn/sXeJm7WMHXzRf8nTtWNwpTxuYHaRDWPiH9XTNK5q1\nWjuduzFaoIihhLNJc0VzO4coNwlWK4EihlJ5r0BwiuYR3NvJVLNgdJlcE7b8LUmA4oJu0bwfYOCh\nT1eOg/L66VePS2KiIWd8fCoNCwTWnKIZX43SSNc/BjjYgFWHjX8viGdZDB53OrhoG7kmGllyXrsA\ncLbq0LHv6ppMeplOOLtoJkVYhqi6oQaAtdU8cZK1Yk2mX2gylBTNuStozxL/R9gEHf/9Y4xMcu2x\nZOqpEkUVUjQfnTfZIQo34C3Ai/KzLCQam2/I+I4aRVnRjGcTWqnnanCVLBPHTFlINNYrpdOsHizT\norlwB2ilnm3Cw3hemn1WGhINkoEXm+59n3iO3GaHKFEWGzHuyqVP8XBv9TDvtmGn7CU/GLpqz2iK\n97xgBfmEfjGkA19FoV6pJuw2KPc9NUem+X8ZPLxxupjm42E6nLpe6/fMtixDuqKRrlb/7zjxz6Na\nvMfbWFTfv7YQ+swK+0DvT6qb4d38Hk4lZVKl8IiUFc4lZM5F6Z+JIpNw+vUo9YtRZBJOr1P1fx3w\nfvV4rozgAAAAAElFTkSuQmCC\n",
"prompt_number": 3,
"text": [
" \u239b 2 \u239e\n",
" \u239c -(-\u03bc + x) \u239f\n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u239f\n",
" \u239c 2 \u239f\n",
" \u239c ___ 2\u22c5\u03c3 \u239f\n",
" \u239c\u2572\u2571 2 \u22c5\u212f \u239f\n",
" \u239c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u239f\n",
" \u239c ___ \u239f\n",
" \u239d 2\u22c5\u2572\u2571 \u03c0 \u22c5\u03c3 \u23a0\n",
"log "
]
}
],
"prompt_number": 3
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll = simplify(logcombine(ll))\n",
"ll"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\log{\\left (\\frac{1}{\\sigma} \\right )} - \\frac{1}{2} \\log{\\left (\\pi \\right )} - \\frac{1}{2} \\log{\\left (2 \\right )} - \\frac{\\left(\\mu - x\\right)^{2}}{2 \\sigma^{2}}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAA3BAMAAADzpFDIAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMqvNiRDduyJ2RFSZ\n72bxr6VbAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAHr0lEQVRoBdVYW2wUVRj+t7uz992umBATwC4B\niZpoawQf4IGV9MGaRjZAUB4IKxIu4qU2KhiIzJN4QbqSJtKnLmJ47T6VRDQd8UI0XCpeYlDDGh+U\nxMRCKCBU1//cZ86c2W5fVnqSzvn+7//+c87fOTPznwWYza13b3E2Lx8gYrd1ze4MonboxuzOIFOJ\n35zdGQCkJgMyiFSEw8oLdFv20XLAsnoVf1TB2xC5FupZXSSvzOz/+bSHamohGsrmkMg6GivMT4oU\nLaPXlwXr6cNM4uGmM9h406nc/v1uQ8ML0f4QjmgsM60LtN/PnpKNVYMovm/mGfDxDKMFUo1yDgOE\nVpx4xRgbcwj95EmWQWrcL7L2ds44AzGef7RAJlEIdAGgM1OvTxgVY2x1CZYBbDGJhmecAYjxTMOZ\nuYxt5ilrjQc7n2MuMeOYbZC2JINjhokVtUJBDWX5l1pkEHU0ATFbkgF5V4bGIZEnM6q24dDbg2WA\noDcpQLSPiUUGWb6b1BCIWpLB6zhRpARRB3vVrEJ2R7QL4DdFaaijwgiRgfW3JiBmSzLYiROFHdhU\n9SwgUoxMpmyAUQ/rNs5wQ2QAD7u9HLckA/JWx+WPFbFfu5W0hxBZEKa7pIMvxd+d55TMYNj2i1qW\nAS7/M3369jJhAjMI/csDZAYDNMA7TEsyILvoS4Atce/c0JEDvC2Bb6qEODTIDDL0pnlHaUkG5Ene\nBvFbEc/cqeIwJDCDwCc5Ns71MgPJuAZqSQarcNNfg/SNea6JAUaq++AgMqc9rMsQ//H0A9ceZHRE\n3BWlOnlxc01ZLnQ3x7+6OAbVeJqLF6C02PS6cJ+kpnrK5OXvamuGNjxRQ3s55+YsdTkJ3ORoBGRv\n6oywfcFxseHackKj9aFFBY3BAvTAqZcAFmo8fplsCD7Kx8VesTZrgSPejNEbMn0QaJQvWBXdXdqw\n0lyfl5CBZYAHwdUVwGJTa1j+Z3yrEZqI9PQLivdjVY0Aa0pnpK0HD0qPQpJiIJb3ElhjthUhg8VD\nwetAazcM5HwkJ9ZKh76IYX/MFSnWgR78oxRslEgDegZYgGZK0HYZDMXm/uBXJmYnmr6Ic7bwyP66\nRDrQgpMlKQj7biV36Rng4xqeJBmAv9iM1+R4OsjaksFFrBt6HOCD3kvvE9KQwbdSrAMtOOLgyaOO\nzQECjQ0zsHoGK5BcfHaoDxWsxIxi+RhcbBoHkmQ/hLZDrJychBcShDScvn6QYh1owbEaQPe8w5fe\nxEN4l67lNmYw14bXYKQYrZWRI58tfMkXGnyiqCL40g/RAoQut5VgCRV1+qXP+ynOaMHRKoSqEfs9\n9KZKAUGYAZ76OnJLAKfEtpPqtuE1uNikksBLP3RUAa6mx70ZkL0gjqEsA1UdgvRpwZkKThODe/Aa\nmsCL9QgpKJ/JIWynMcSdD91Cs3CGZ0AKUEj34YWWavosNKzRBaAfxioA1+OXrZ/IWHCuSDv3Zavb\n8GAtmGYwB0qooRl4tNyI5VP4hWzvG8ixX0FoBoeIM7DYNA3j4ti/cSo+OEj+VQCrbdq5L42eA3ID\nZTDuIoCzQA51rteSeyh5D5y3eh6lDrKLUg7go3PMK2za6oeMA8mbqQqPMHwPGr2LPMGxGu6cScD/\nccMn+Wd8cisf8fnIkzwf4HDTT3Lk1Bc8lHe7IPQNhAvJ5XfYlPmrwh2qC/4eaMGRPL5NS/AdhqYd\nFe9BsT5YkIOlcPTIXZRfBRDf0n2y1KDY9MTDY7DIQ8yp3wsHuu8HWFn/xyGesQq5ultwVaEHk60T\nq8FijA5X8GJooc4rdvzzVVgG1et4L+jeQVgvqWLTEOWmtkOH7bYFjlWtd3YQo6nKjhaTIhRXLYK3\nS+5diTiIn97FHjRu74GnRsuIo7YQjHNg6lwTXoCBiklyHMnficNfXafw7edtrJhUnAz+U3I9EnGw\nBiLk+RAt1YfPCv6p359VsSk0qvdOOFZUHoXW5wD2EDPTpUiG/CccVkwqnQyO2Zy0fKPgR+IXFQHW\n0+TWEUKUaKrYdMk49E5I94pPZHUP9WISOOyE7mvzMayYVDoZLE84MhUpwldyJ52AM+t6Bj+lEE84\ntIlMuOnpPBMGllw8RB6P5RDi3CkJUUwqQqKDHP0hGQFWFr0ZCB54AZq1JeMHvHplDrVV/ULCkG8+\nadaZen2KooEy7bwXUkzOtL1YnGmES48TsiUlC0n3vXRJBDzPwNyvR3fbFA6zjtHiOlIQqOk+ebVp\nqUGIE7Il3Xni+DT/iXM0HKuk7AQb6FXDePg7zoyb+El5xoE0YBvwJV2s16cZYaRCBOkuwCMDbaaD\nPi0mmbvp66mmlQYhTuhakkHgoqIOMdprkGD3wPjrOy0mia75lsg3r/UrcUK1JL/bw7AVD+QgnKc8\n7zwaVkx6qGmNjyFuTysKEpAJ1ZKCVIKn5cqmIixgD8yILRyqn0+LSWU3gUIOJOwmdGYJmVAtyaxR\nLP1ox6rxZxm1V3kE4sWkMJvq3+juXdKU0CSiE6olmSRuLuygZX01yG5B6nu3j2FeTPodDZiV9fq1\nBu7GLjqhXFJjLXrj+NuNahsrCs8aNOpeqWETud23J07n1bpoBazM2YLuUws9quBsQumKWK2FZf3t\n3f4DvSUcrYcgV58AAAAASUVORK5CYII=\n",
"prompt_number": 4,
"text": [
" 2\n",
" \u239b1\u239e log(\u03c0) log(2) (\u03bc - x) \n",
"log\u239c\u2500\u239f - \u2500\u2500\u2500\u2500\u2500\u2500 - \u2500\u2500\u2500\u2500\u2500\u2500 - \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" \u239d\u03c3\u23a0 2 2 2 \n",
" 2\u22c5\u03c3 "
]
}
],
"prompt_number": 4
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's pull out the derivatives of the log-likelihood: $\\frac{\\partial\\mathcal(l)}{\\partial\\mu}$, $\\frac{\\partial\\mathcal(l)}{\\partial\\sigma}$"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll.diff(mu)\n"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$- \\frac{1}{2 \\sigma^{2}} \\left(2 \\mu - 2 x\\right)$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAIwAAAAqBAMAAACTq67+AAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiK7mat272aJ\nRFQidGHIAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACXElEQVRIDe1WQWsTURD+NumabHaTLj2IPTUF\nr0pADBTaEvRQ8GD2Yj2I1IM9eTAq1CoIQTzpJZdCCwXzB2qDF697ERUFC/6ABAoetaFKb67z3m6y\n7z3fVpL14KFz2J35vm8ms7NvlgCqzarAWPFyf6w0Jen0+j8pA+ukjDLZODyZTTwL1fuvZpN7fvhC\naHCiIwRJrl1LYgb4MnB97fYgku4C/kwiNMEsjCY+72kYEc9oeBGyfORdFNoiFvkiTrJjreChUEH+\nh0Yk4kZbIxCgL0Cmry8j4eeFHI37jmPZPibasJqqgPDpzZvbLSDUqfwwvsO9rg+ngmxjCEdO1zd8\ney5bBl6rlBwv8rBKj9bATE/mgCoc1+mXPGAfUxeYnQUCxQ4oa4dl5m6B1fjukj9Uh7iBDJFAnV2S\njZfZJJ5q3FdlDMdki13/UoY9VKmBa3gAzJssITaOo14DdUlvVGfO2j0OsxEvARs0BvOXIysZXnJf\nwqIyCSN+i6c85z1gzn9arxiHyP28KpXheLd3FzcIfihRw+Ai6h4L2PGjyVdKRx9a7HwIxvEzW9Mf\nOwSeEwhhZXewsseYLC8G5MuC7k/XbMeYuLL81RBlN0K+0Ip1Gs8RaHFlgblQvRreVmqa5Biail15\nlZ2oDfpsMUt4oSEJRD/GQ2ll30QKszOQHnO3PZmklTUeB8ERin6xJlOjRLTKl5/sr3q49O2rO0qi\nrK3CeAX7ANgNApkZJaJVzpVRTPs/iVZ2sgOLukljbGXpgGSaaYqEqzzj4kqK4VIDfGVP9cyFdM3w\nlTUebY/VzG/DEqMiKJVV6wAAAABJRU5ErkJggg==\n",
"prompt_number": 5,
"text": [
"-(2\u22c5\u03bc - 2\u22c5x) \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 2 \n",
" 2\u22c5\u03c3 "
]
}
],
"prompt_number": 5
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll.diff(sigma)\n"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$- \\frac{1}{\\sigma} + \\frac{1}{\\sigma^{3}} \\left(\\mu - x\\right)^{2}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAKMAAAAqBAMAAAAt/fajAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiK7mat272ZE\nVIkJo9jjAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAC1ElEQVRIDZ2WO4gTURSG/0kyGWMm2WDhYgo3\nZsVHF1bWRbAYtFsLAz4KLZLGUogWoiK6CBZqs9uZJqJgISgGVCwsTGshRkQ72SwibCqzuCAiSzw3\nk8nex5xs4i0m53zn/P/cO/dOEkAae6RYDjku94THZ9fG4+HdMt15NdyS47KWixPhluA45yNxTspx\nScqFnJTjnI/EOSnHJSkXclKOh/ikr1xSKCfluCL2kx+4rVBOynFF7CezOFmUMSfluKztx89QasiY\nk3Jc1g7i5cwgpICTclzWDuIjg4gC596v+3IexAo/EVDz85VA6YpZGE7sMl/f5lHtDV9nKhHlSWlN\nOcBtuZ5Gg9RSjkJAgdpmaEZUPPbtK3fT1KKpEORJOPZpFnjR7XIdjKVb4ASCR5rDqoxlbz/ngQeq\n1qreed/earcZy3idvL4Dc6rlbrxcnAGSORWrGWMZFUvbD2tD7V7Bz8wXINV/LDsOibEX6PZGx28O\nLCd8Kq6isr0B4ZcSvyjWtBDOeBRmsI+uiPXFItZGOp+ffpzPVzRMqbAkPyenl/4IMMSSqsEsNalY\nuFNAtFVUC7F1WHTOCypVM92yvZoRDWJ7omWUmi1KNsenRAdLY26PU3fKwiG9AEwVcOt0U2TBcDfi\nHZtu5lQCEvapzdJpJTuiTSxt+WY1QYdQGtbdpdpbyiMNCRqhZkkzyPV6ZqH9uMjKc3JixC6tQxlZ\nr5fSd9dThcuJOne5EhbvOujTeNFeD6sLZuXoMsbwtwd2mT8p8eIYftQa67+D53nZa7NkXe92fxt4\nkr4CK1nP/mtURgDHb3y+YMw9Nvfh8BlEi8lHIzjoLdZz+IdPKWSbeEfLbq8aN1PawhM6ea75n/Uy\n8DETLtiaTtRBr6o+6BE+1NnIeclDZEHvFqfwgA5HzqcymDeWaK0h9j8b49813rSPmve/iFMNk45I\nrGs1Y5LAZHVlRP2g7R/qqadnvG7MtgAAAABJRU5ErkJggg==\n",
"prompt_number": 6,
"text": [
" 2\n",
" 1 (\u03bc - x) \n",
"- \u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" \u03c3 3 \n",
" \u03c3 "
]
}
],
"prompt_number": 6
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Second derivatives:"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll.diff(mu, mu)\n"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$- \\frac{1}{\\sigma^{2}}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAACgAAAAqBAMAAAAzGARVAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiK7mat272ZE\nVIkJo9jjAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAArElEQVQoFWNgQAFMAihcMIc1H1OQMb0fU5CB\nYf5wE6xYb7EBM0QGr4iQMQioMDD8RwAKXPs8HTNC+QI4FDCMZA5g+YohyHWA9Q+GIAMD70csgswT\nsAg+xyLGXYAiyFj1//93hnMMR5BFnWsupwaw6N5NRBJkXMvA/YGB6///D0iC7AoMPBiu4d/AwIms\nCKw+3oGBqQFJJ5gpL8DgJYAuyHaB1QZdjIGxch5MIQCENjK5UVh8EgAAAABJRU5ErkJggg==\n",
"prompt_number": 7,
"text": [
"-1 \n",
"\u2500\u2500\u2500\n",
" 2\n",
" \u03c3 "
]
}
],
"prompt_number": 7
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll.diff(mu, sigma)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\frac{1}{\\sigma^{3}} \\left(2 \\mu - 2 x\\right)$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAHAAAAAqBAMAAACZ28LpAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAzRAiu5mrdu/dZjJE\nVIl1eHuDAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACH0lEQVRIDe2Vv2sTYRjHv7nzcrmLuZYuDlZ4\nC2K69cC6iGChmw6mUKuCQ0bRJSDUrQZnhQwBEbSc/0EWQQepq1OzdjJLE3DQEulQpcTned+7y+W9\nS+A69xme93k/3+d7P5434QCOBZnzp9Vhfg87rm6d0Qj73Dh95OfDmT4b68nJ01gtB3E5vXA7aW0V\nuLF5N82JJPizdMMCjDp+f0sLE7yQ0u0eigKlZkrABKc2LUo+Si0U/2mYt0luNPWG70BhmG2c4Fd0\n46EEzhDlJuy6rhKvbt/aaQCqL6HfkfV+D14LTi0hhNzouW+deeCjLn2QoE0PXMNKV1fb8IQ3rPjA\nAUbj4LbXnKw1sOuHoHrpHcdlqiQ3UFjjco9TMqRxmwi5NpIC18wx1+CcMvKjVmq4jkfArskt45Ac\nex0IgKavYtCnHQUP5xPwHG2Yp55EcWJeES9hi/FwrIDfi6IPmLu/tlrGCaw/XyWKkuT73Xu4SeR+\nSK2eeyRL/gHQuFqV40GDzysRkl97Uf0ZEFyMBWtelo6vSFFtY1krzGYMljuydGuKlBqxklV4sVy9\nFOrrar2trpNlYrY0FsLhgP7IHPG41VbP4eUZl0+VaAZqnZldn2Tj4Wh0vNwx/85szRC/PD5Y9x3f\nfZWhzULGG9ARlgd9vnueoPO7eKbv6lwA+yjPnaJeOrBCPdrkWVcEPos8hqj3Qtd8H9W5VuPBjshl\nCJv/A/D1jnHTPoltAAAAAElFTkSuQmCC\n",
"prompt_number": 8,
"text": [
"2\u22c5(\u03bc - x)\n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 3 \n",
" \u03c3 "
]
}
],
"prompt_number": 8
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"ll.diff(sigma, sigma)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\frac{1}{\\sigma^{2}} \\left(1 - \\frac{3}{\\sigma^{2}} \\left(\\mu - x\\right)^{2}\\right)$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAMQAAAAyBAMAAADxQ2lcAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAzRAiu5mrdu/dZjJE\nVIl1eHuDAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAEfklEQVRYCZ1XTWwbRRR+XttJ7LWNQYCEENIg\ngZQTWalF4ucQC5QDB4hBpfSASIBDEZdGPZBbHP6pKKqrUkokBJYQB8qh5lY40O0JJAS2gAM9VI4K\nFeJCQ2lBpUTmzezMzhvvjHedOex873vfzNudnXlvF8DSaqGF3AHlLTsHXXR6JnSc1foc0xigtkwt\ngvN//NonZjos15Umv8gU5P3PhkU890L1X2JmgI9JjbdwgBG5d5gYBmz14ZpBpBqz8VO3GBEXmsQw\n4Deh96dBpBrVdSUxQgyYopP9pAsFH6o5jBDHFWvpc84ntIg5NQikg4Yo/+VQI/3VuPDmsNqeJzhR\nUvdEQ5Qaptawam8b5hjjF3iWe8uXpYaGOBNK0tq9xax0kjwGZ/gaef9IFw2xLylXzC6Axb4yUvrD\nMB9yyY1SR0O8ITlLN2Sw2LHwdmrAON+S75uE8Mec39cAPpIj7LMa7Alhza9EJAkxPWZDnYSqWlpj\nMqtRawp6piG6vUfe6wiAl8I6Xhyt8uSeMHZ9GaME+Iwzn0R0YT3qyVVGhZsIZ4P5ho2NuKkuQGWz\ngldM28lVmYueb7faz66JcszlQb4O8Pn5n4Si/F9C11vh1F1700KsJUYSAp1HhkNB+Mm3N+gLz3Ra\niFfJjAk4qxnvqsYStboCpIWotBMjCZGL7lMwfxM+gkuB6NNCiC15GuAlcwJv9fmvf8N33NT0FQ0l\nyhii0EH97wAjqfce+HTlKCa/up73ZQ0lkjU27SlKfCluAW/bnOACXGI/AlTbmj4EQ90EfSByRiG8\n2zawHe0id4PWbaE5E2IW3Yaq2BV3ctXG7UgzuBmv4HOJbO8rEPdGiJhNAB4C5y/WRz1ij44PscTE\noCwLVWxDaTMwY/jXwcOj3dbshoKyDgJcisY4Q1xcYHwMf91YH+f7m9yK2w/TW3DQfN2HlFPWQczv\nXUG5QtSCqToX8C/GuTYsPdgXcnmpbBe28h2AYlOz8Y6SdRDgVMi9xReuvahVBJUCX+Q1vhSDZ1an\n8RCQ5u0/uHYS7VyoyfhcqDoIA+LVuhjNhPkorx2DqP7HHgru0wZNIAMm+N6KdltRtE95OXjT6uck\neTaaBk9EA2Qyd46GUnQPhSB/3SXy6tpT1YdTpZUZ4tZCguTPR75BdybxIywE2iYlSdZBKGxptw2V\nm5J9wOaNuI+Jawrn8/YNh1fjOgiu3QpwN1aZJnwP35IJ0iEv1F88fe7hIK6D4Ou1M8f7x797dxf4\nd5x/yORTLPzI8d6BMj6LqoMArk81/Bs5helvOExZyZGQ+KmGiaxiFNGlEY0yHwfoMWVk7x/BJN0B\nzCuk9UJiEIhfiS1iZoWYeee7kFum+viPgJIA/BTcalJZLP7xP8fgNKPi+I+AkrjxLoP/ikllsfjt\nF/r5D0zt66aprEfh/lDh7H0vwLt7ao2ZI2SyMkk8FqsXRqkM9oJVM/Evo3WWiKza19ad3sbM5XDN\nhnbHOTu9E9a+Tljylncym21MuWFjObff5ZiUP+scUAydrokc3nNc/j+rEw/q579VVgAAAABJRU5E\nrkJggg==\n",
"prompt_number": 9,
"text": [
" 2\n",
" 3\u22c5(\u03bc - x) \n",
"1 - \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 2 \n",
" \u03c3 \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 2 \n",
" \u03c3 "
]
}
],
"prompt_number": 9
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"simplify(_)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\frac{1}{\\sigma^{4}} \\left(\\sigma^{2} - 3 \\left(\\mu - x\\right)^{2}\\right)$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAALoAAAAqBAMAAAD2RlyBAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAzRAiu5mrdu/dZjJE\nVIl1eHuDAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADw0lEQVRIDb1VTWgTQRR+2W2bn022RVBBLWwP\ntfZiA/biD7ZailAsjWL9AQ/x7MHowd7aoOBBKy20WgSVFBXEi/GmB2k8CIpog3jpoTQq1epBSkoV\nf9r43uzO7EzSBpJYH+y8733vvW8mM5NdADKPhcNamd5jrZU0gNZ5Zg3VAQb+o/qDf7FNWlSoKGv3\n1gm+EvBMNCvqLyzBVwK8Ud4tq2tXOFu+n+20ADbxflm9JsLZsr0ZrsHdbcMZmMnqE5wE88hhJ1+i\n84eNBYCaqNMmq68XUh/hlMAlAV9K/wNgZAvVdZzVsRGYCHNcog+R8kanSVp79bAQGoS2lAhKA/4Y\n1t+3e44OXU/wbh/x3GYsjlRvJNRYigJJCmZpWGFtbSxLSbQx24lx77duhlsFUwhOIxWIEO/L0KjY\nQykyWZFL6BGYSVC4zuUKkAeZN/AKR7qWeXZCiu9JmKB3CWqHyWdwWM28GTDqp/dj2jtfUDPqMsFM\nMOlGiIKD0ELr8YUVWg20OPhyORIOudePlyxyAPBo+p3lRjZiO/Myn1Xieh4ZSxwJ/4OhxqFcLoKP\noDlg/zW6EUYc3PeUnW06v6c/5twXYrRfNu+ONmNcez26wyUFMr5ECHfhYw6DnwUiqWUCY/46gHHB\n2CsVIYDO5tuWBvnuSPm5JAa38PFEoCUtJXA+y8yGwgBTgl2GnGvE2nt1CGDS4kUNN8g226HvAvpB\nfFCa/dfcrAaediqaoIGZdIY2YbC14zDgVMhOD0MVnRSpozSuQbXaGMVF1NnO6L8BNqiNLKrNQhXN\nTjuD/4vbel7NRBIsAPdGLbv5dgY1OgktCwbtQL752sFPd5hOdRT0JVMpCFkD4LXkU2X7wGrYOxMR\nm68bdqYYqw6BBMzQj/+EC/gJ1QtPlfRkugd2IdPLWZ120bbtzrv+EoWNfR8cWnXPew8Sgb89tPg5\nRpdbsq19TV8TGG/hnIGfEMfOOer7OFHE+8MrvaGcBj3OO/l24E1866jP8VwRjy9Y5TOglJriF5lx\nnjDxJJmxS8zZ1XwHKJ8BpaxBRFVRDu9w9ZYEp4r4Vuna5Zd1CKI24UAtxtXdCUVVIdB5X2EqEBYc\nffO147ncYqi3d9HuCMRFtlLwHgWenJzqwEl0Z98195JWqn4Zl34VAvOoc/d72lY7W6ko7w/FAarr\nIOgs26Z9SZ6u0PsTAHiw6tc1UFehKm8fR4AX1xPlBPMHlKjsQL+IrS0WPLYUieZ/szXNdI5Vaf2m\nIg76iBqXGXVRn3as38rr350XlxVqKWz7CxhC3OWXjyL0AAAAAElFTkSuQmCC\n",
"prompt_number": 10,
"text": [
" 2 2\n",
"\u03c3 - 3\u22c5(\u03bc - x) \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 4 \n",
" \u03c3 "
]
}
],
"prompt_number": 10
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's find the point where the gradient is zero for $\\mu$ and $\\sigma$"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"solve(f.diff(mu), mu)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\begin{bmatrix}x\\end{bmatrix}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAABYAAAAZBAMAAADHzMPdAAAAKlBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADmU0mKAAAADXRSTlMAu90iEHarRIlmzVQyXQFQhAAA\nAAlwSFlzAAAOxAAADsQBlSsOGwAAAF5JREFUGBljEFIyYAABdiUNBgUwC0SYUZPNEppc2go104Nh\nUoAWlN3C0OuwGcp2YFiIZO9VBJvnDgMLVM0W1gsMIRA2+23GC8wTIGyW1JCImdR2M8gBQGDGoAQP\nHzUAIM4WtLB8TSEAAAAASUVORK5CYII=\n",
"prompt_number": 11,
"text": [
"[x]"
]
}
],
"prompt_number": 11
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"solve(f.diff(sigma), sigma)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$\\begin{bmatrix}- \\mu + x, & \\mu - x\\end{bmatrix}$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAJgAAAAZBAMAAAA/APpvAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAu90iEM0y73ariZlU\nRGb40ErCAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABcklEQVQ4EeWVMUvDUBDH/6YPKW0j4icIgQ7q\n4iQ42cXFpVqIi4tbP4Hg2KFDdRXBsYiDm84uFnR3cHEsfgDpUIpDB5OXu7y7QJYQJ9+Qd+/d7365\n5AWCjXAHlQw/3EJQiclK9v63zJQ5laJ31hqUOJc/lPWBy7Qj7syMrj8fck06ihNE6c7ugG0tO8TX\noM01NDuKE0Rp2SnMUsvu8dh55RqaHcUJoqzsIEzGSWJqzbWsg2Ou4Jkp00uq2mfxPlGqs9jUPIpz\njSjqnUfRky3/sVdxYUpsAZZSsmYAb5x+YHwAWF3AqDIIKkuklJJ569ifjC3Bsrf6DMOsxgaC4gRR\nStYN8LI7sQTJ/OXKrDZFN3l4HoKiLaKgZO/Pozp9VCQzt8Orb8CjY7HFgiIZUVp2wzcH+DFp58Nl\nICixm4SqswuX9KcujiO5EpRitKy2yCWzpRGyYkp15gdZdS5YE+tiSslERblQvbNyCldVsSys7r+5\n+QtCEFzfw9ZCtQAAAABJRU5ErkJggg==\n",
"prompt_number": 12,
"text": [
"[-\u03bc + x, \u03bc - x]"
]
}
],
"prompt_number": 12
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"What if we want to use these calculations in numerical optimization code?"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"d2f_sigma2 = lambdify((mu, sigma, x), f.diff(sigma, sigma))"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 13
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"d2f_sigma2(1.2, 2, 3)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$-0.0463620287292$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAKUAAAAPBAMAAAB3ghJhAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMpmJdlQiZu+7\nq0TEZSulAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACyUlEQVQ4Ea2Uz0tUURTHP29+6Zsf+rBNIKEL\naxGUA0EQGM6mtg6RbtQaKkeowEGyCJEmqE2bmYqCdo8WSUiMVgslkhdEBFkZFQQlDkG0tNGKtHQ6\n997n9A90Fvd73v2c+51z77tvwA+rZySlUl9DaYgMnoCb2ecymz2b4kbfUx/b3XMprHMDZUL7e7LZ\nGimiibLREchb+1Tia9DD2s0h7DQFl6ZyuNWaZ6hs8DYiG4ThA9FqtbpuSLNrz6CJMZRxAo6pB19H\nPQIZ5qj7RWOSHcS9hENs0uA38IgrcI1wEeYN+Qif0KTm+RIKjjwZtW97lKSehllaWqOrksWSJL4b\n/AAWU6dhCxE5oZQhyzDkaFLzXIelvDwZjQQ8Xvis4MYnJQ1UlKfG7Y54Lj3hpKr4aoj1By6UNfHX\nYf0Qz5y8A6N3xXP54EhZ4e00tg0eUVmwsllGp5Oo3nHVpKeGYIX70qcraafarQ77J3RlwKiVC3jW\nsss0hI6nadlJLC9lS5nNsgYpP7qsVkdcGYQwBO3SlCJ+bBYbjSCeVYdetepUqqVCeEqyDv8npS+P\n6OH2WZlr0QYdcja5etWnEJr2qGizxF7vXesl8WQFSmVZEZuSNx/dgLDcL7+sD+4R/S0/OSwVisD4\nmaE8CKmFnGNBplB6MaM85YBKrl0k/ivWSnQNxnycp24epPEuWSAdGiKy6GiiZnQIK6mdKt26sLD4\nOCfvvVRurBBfk/cufdanuaWxw2Xsq98g4WLJRTBEdBeKFLWfGkaRr+afxjx1XXsRDa7Wy3lOcgDO\nmzL5cOuKMxAvEpJ9GTJRDK2iibLREcxbD+V6oxUaPeIZa5qoSyHHM5rz9t5sf9Lg69mBt3x2GJcW\nxdOQYafJRZNNS6zBsRTBeaNEXq+7DPSn4Ev3K7kw3e8IyKedNLi9Wl0h9F7+Q7DluhkSUYWa1Dz/\nY/IXed4JxwCefm8AAAAASUVORK5CYII=\n",
"prompt_number": 14,
"text": [
"-0.0463620287292"
]
}
],
"prompt_number": 14
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"timeit d2f_sigma2(1.2, 2, 3)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"1000000 loops, best of 3: 1.61 \u00b5s per loop\n"
]
}
],
"prompt_number": 15
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"from sympy.utilities import codegen"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 17
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"codegen.ccode(f)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"metadata": {},
"output_type": "pyout",
"prompt_number": 18,
"text": [
"'(1.0L/2.0L)*sqrt(2)*exp(-1.0L/2.0L*pow(-mu + x, 2)/pow(sigma, 2))/(sqrt(M_PI)*sigma)'"
]
}
],
"prompt_number": 18
},
{
"cell_type": "code",
"collapsed": false,
"input": [
" print codegen.codegen(('d2f_dsigma2', f), 'C', 'd2f_dsigma2')[0][1]"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"/******************************************************************************\n",
" * Code generated with sympy 0.7.4.1 *\n",
" * *\n",
" * See http://www.sympy.org/ for more information. *\n",
" * *\n",
" * This file is part of 'project' *\n",
" ******************************************************************************/\n",
"#include \"d2f_dsigma2.h\"\n",
"#include <math.h>\n",
"\n",
"double d2f_dsigma2(double mu, double sigma, double x) {\n",
"\n",
" return (1.0L/2.0L)*sqrt(2)*exp(-1.0L/2.0L*pow(-mu + x, 2)/pow(sigma, 2))/(sqrt(M_PI)*sigma);\n",
"\n",
"}\n",
"\n"
]
}
],
"prompt_number": 19
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Too lazy to compile a C extension. Do it for me!"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"from sympy.utilities import autowrap"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 20
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"fw = autowrap.autowrap(f.diff(sigma, sigma))"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 21
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"timeit fw(1.2, 2, 3)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"1000000 loops, best of 3: 214 ns per loop\n"
]
}
],
"prompt_number": 22
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"logpdf = autowrap.autowrap(ll)"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 24
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"logpdf(2.3, 6., 1)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$-2.73417022465$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAJEAAAAPBAMAAAAWmjAjAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiKZu6uJRO92\nVGZ6zyUAAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACV0lEQVQ4Ea2STWgTURSFv0mTpmmbdHClKCSL\ngrixoqjddUAQXDUgdiNIEEyNmwYEf1aJP5QuxBYR60DR7BWsWumiiFFQFBedlbhr90KrFrX+0PHe\nO5Om7r0w5917z5kz8959EEeqfNHT9Mi470+l/OlAi6PyJFySgxXfx6m88TBZS6u0/84j8TJ5TuVR\nbCfzS7PVMAzrI2Q3tDgMqZpLTnrrJALnICZraYXe1ugs0BOGTZVHcQUWNDsDXSwGrEneO49zYdKl\nswh1PsEtTBZrhWYP6SqJu68jE8OnMOlJNgW3edBwvkp+arfAoksGMh6XYMk1WawVOvdNFPLprTHg\nRk7Q3ZS+7W4qdpL6IazDo8BksVbo9Kx6/OskjSFX26QVEiXIFdtOVZzv4tQUymQCSvf1z8hZd63I\nLLZE98+ouCfLsUMCu9h0yiyTEnpY9m4yBaXze+kJSLspm1bLrKNqmfNZl8ycHtmmU16GGDuZTEHp\n/Bc6H6v+vIzxgEa/FKPakYEVbFlws17b6RU44qS7M5mA0X37yNnvjBXtLYPeepTr505ALdhB20lu\nDnJOSwEmUzC6p0Bug6twvxG9rfiBlPnm6xC61JY/Tkz8vmG3AEfmptd01Y1kqjVaZif/9ATGXPWw\nSJboNafhEsiVG9JCVr1PJOV/eA+nMVlLO09WzmkWeWO/COI47t+5Rr4gOyjBCtk/2v8hjzpl1akj\ncOYwmUFEX2ZnwAi5Z1LGMRCGa3TU4WQgky6PNqQ/Hd7k7ItBuQHPpXJm3nqYzCCiM+XrMtbKuHzu\nv8VfhqXApTA0Rd8AAAAASUVORK5CYII=\n",
"prompt_number": 27,
"text": [
"-2.73417022465"
]
}
],
"prompt_number": 27
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"from scipy.stats.distributions import norm\n",
"s_logpdf = norm.logpdf"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 40
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"s_logpdf(1, 2.3, 6.)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"latex": [
"$$-2.73417022465$$"
],
"metadata": {},
"output_type": "pyout",
"png": "iVBORw0KGgoAAAANSUhEUgAAAJEAAAAPBAMAAAAWmjAjAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiKZu6uJRO92\nVGZ6zyUAAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACV0lEQVQ4Ea2STWgTURSFv0mTpmmbdHClKCSL\ngrixoqjddUAQXDUgdiNIEEyNmwYEf1aJP5QuxBYR60DR7BWsWumiiFFQFBedlbhr90KrFrX+0PHe\nO5Om7r0w5917z5kz8959EEeqfNHT9Mi470+l/OlAi6PyJFySgxXfx6m88TBZS6u0/84j8TJ5TuVR\nbCfzS7PVMAzrI2Q3tDgMqZpLTnrrJALnICZraYXe1ugs0BOGTZVHcQUWNDsDXSwGrEneO49zYdKl\nswh1PsEtTBZrhWYP6SqJu68jE8OnMOlJNgW3edBwvkp+arfAoksGMh6XYMk1WawVOvdNFPLprTHg\nRk7Q3ZS+7W4qdpL6IazDo8BksVbo9Kx6/OskjSFX26QVEiXIFdtOVZzv4tQUymQCSvf1z8hZd63I\nLLZE98+ouCfLsUMCu9h0yiyTEnpY9m4yBaXze+kJSLspm1bLrKNqmfNZl8ycHtmmU16GGDuZTEHp\n/Bc6H6v+vIzxgEa/FKPakYEVbFlws17b6RU44qS7M5mA0X37yNnvjBXtLYPeepTr505ALdhB20lu\nDnJOSwEmUzC6p0Bug6twvxG9rfiBlPnm6xC61JY/Tkz8vmG3AEfmptd01Y1kqjVaZif/9ATGXPWw\nSJboNafhEsiVG9JCVr1PJOV/eA+nMVlLO09WzmkWeWO/COI47t+5Rr4gOyjBCtk/2v8hjzpl1akj\ncOYwmUFEX2ZnwAi5Z1LGMRCGa3TU4WQgky6PNqQ/Hd7k7ItBuQHPpXJm3nqYzCCiM+XrMtbKuHzu\nv8VfhqXApTA0Rd8AAAAASUVORK5CYII=\n",
"prompt_number": 41,
"text": [
"-2.73417022465"
]
}
],
"prompt_number": 41
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"%timeit logpdf(2.3, 6., 1)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"1000000 loops, best of 3: 197 ns per loop\n"
]
}
],
"prompt_number": 42
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"%timeit norm.logpdf(2.3, 6., 1)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"10000 loops, best of 3: 86 \u00b5s per loop\n"
]
}
],
"prompt_number": 43
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"x = np.arange(100)\n",
"vlogpdf = np.vectorize(logpdf)"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 45
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"timeit vlogpdf(2.3, 6, x)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"10000 loops, best of 3: 51.8 \u00b5s per loop\n"
]
}
],
"prompt_number": 47
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"timeit s_logpdf(2.3, 6., x)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"10000 loops, best of 3: 108 \u00b5s per loop\n"
]
}
],
"prompt_number": 48
}
],
"metadata": {}
}
]
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment