Created
August 19, 2014 14:13
-
-
Save wiso/bbb0f470668e3ae30859 to your computer and use it in GitHub Desktop.
truncated rms bias correction
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
{ | |
"metadata": { | |
"name": "fraction rms" | |
}, | |
"nbformat": 3, | |
"nbformat_minor": 0, | |
"worksheets": [ | |
{ | |
"cells": [ | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "from sympy import *\ninit_printing(use_unicode=True, wrap_line=False, no_global=True)", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [], | |
"prompt_number": 92 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "x = Symbol('x')\ns = Symbol('s', positive=True) # the rms\nf = Symbol('f', positive=True) # the fraction for the truncation\na = sympy.erfinv(f) * sqrt(2) * s # the limit of the truncation\nnormal = 1/(sqrt(2*pi)*s)*exp(-x**2/(2*s**2)) # the normal distribution[0, s]", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [], | |
"prompt_number": 93 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "# check the truncation range\nintegrate(normal, (x, -a, +a))", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"latex": "$$f$$", | |
"output_type": "pyout", | |
"png": "iVBORw0KGgoAAAANSUhEUgAAAAsAAAASBAMAAAB/WzlGAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEJl2IquJVETdZu8y\nu83OyatpAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAVklEQVQIHWNgYGAQMnIAkgxh9QVAkvGHBIjD\n/gFEMrBNAJFc0VMVQDS/AYhk8D8ApurBJEMmhLoKoTaDKcYfYIrjI5hiuwCiSvkWgKi5kiCSYbEK\niAQAn90OqI5y3l0AAAAASUVORK5CYII=\n", | |
"prompt_number": 104, | |
"text": "f" | |
} | |
], | |
"prompt_number": 104 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "# variance in [-a, +a]\ntruncated_var = integrate(normal * x * x, (x, -a, +a) )", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [], | |
"prompt_number": 94 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "# variance\nvar = integrate(1/(sqrt(2*pi)*s)*exp(-x**2/(2*s**2))*x*x, (x, -oo, +oo) )\nvar", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"latex": "$$s^{2}$$", | |
"output_type": "pyout", | |
"png": "iVBORw0KGgoAAAANSUhEUgAAABEAAAAUBAMAAACZjst6AAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAiXZmMs27mSIQ70RU\nq93rZ8ecAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAcklEQVQIHWNgAAHGskowzcAgwXAVynrB4N8A\nYa5m6N8AFWQ4PwHGegRjMBbAWOYwBrcAtwOE7XrmGFTH+v//gUKMSiZQRSEM6RAWxwcGTwiL81c1\nVJJh2/8ECLOHYdIPCCuKgeEDhJXCwKQAYTUp60AYAF0TGH8jlbKxAAAAAElFTkSuQmCC\n", | |
"prompt_number": 106, | |
"text": " 2\ns " | |
} | |
], | |
"prompt_number": 106 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "ratio = simplify(sqrt(truncated_var) / sqrt(var))\nratio", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"latex": "$$\\frac{1}{\\sqrt[4]{\\pi}} \\sqrt{\\sqrt{\\pi} f - 2 e^{- \\operatorname{erfinv}^{2}{\\left (f \\right )}} \\operatorname{erfinv}{\\left (f \\right )}}$$", | |
"output_type": "pyout", | |
"png": "iVBORw0KGgoAAAANSUhEUgAAAUAAAAAwBAMAAACS677sAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAzRAiu5mrdu/dZkQy\niVSnpIUaAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAF40lEQVRYCc1Xe4hUVRj/zb2z857ZSTPMyq7k\nIhHhVYtIAiezCIl2yNRQtEGwspQdEDWlxUeIaFqTqdAf1WyQglBNpZY92OmPzALdgSwpXB2iB0bi\nMrZqWm6/c869d+48xA3buB+c7/t933nM757znccASnyGBTxq9PYmBAe8IGrCtPs7mhDMeWk2VzcS\nTKQ9TjBQ8DjBsJf4ockS7/M6wbleJ+ipTdxkiWMZj89gvPyfETx0xLzqsRo3Sdi46kHtAVJRmS6T\npJ+wo//ONhK8wxngLnXl9TmBKwPtaKdRbZXUKnS0Esu8bCxfrRg8mvXSjvp+q5zezzho0CBs+rKu\nxgkxg34DiHeWcMxVcTVwmN1ZT9po8PZu6O7GY0x6H7L4CsBYd01T/LUVjRSbVqug5tBqKTvNtNkz\nHCzBnY8+VBuwvK9qotoC4d7H0mpwJosEjujzVpQdRwE9JWzbMhOLVaCpjmTt8K02AALHL1UdIi2D\nE5wTl2izj+rX7J3bMfOx9ITNnRnsv6itNdsmltgkB0Ta5wBBZ2TR76AZOePqL6DPoMKS02X4BLiM\n+E27Yr4NgPHvF6oOUcBAaE1NxJcPlD/6oHQQ8bR+CaMR60MJbwyQXCzJhl0siYy7x0JEp7p94k7h\naxdv47eU6qpc7hgbaxzcltM2sGwoh8BfNbFuI5Hq4bwIghU8CYzSsqpBNEU7k0UCFaO+yUEOeF6g\nRB+VtkbA5nLKDsfzNgK+r0KJfBWbYGIahVPx7u5dnD5FsE8QPBk0VKdohnY6i5xKFaMmQcNxJFCH\nZlzOy+jaKst7R1gnP0UKruZJ2A/s3LjQalI14Qq0o8xnSzh7ZOciGPrMqomlCF5gkUz1Od+Z7204\nNXzfxicWXcLn679cno3cvAkdr/MoSrNR5KkdwnzL0iDRsyIkP0AAnofxY79OXio4nBOBWukuYb+p\nl+ygP49sDcH4ZqtKpIo8sCOi8X4jug7rf7yAEWCu4tqsfgahjNwVHILSmhW6S6h6mTySkUSaKlqg\nygGH8U2wSKjX5hsjwDboL6KtKLFQvb23tM/Q29MnN409G9xYRvUU3coBORjkS30LsAHX0xsBjQQf\nBC4gWMFkRsJlKkwpCn0SqP/nBhRPsC7AgmeT9qGQ4n61EhcY86qQGxkAWlIIbFrxpsRXUFyKSJpt\nfCyx/t27n8Z1RIrgvWJ1tAvc8UCoQAW1H0XGNEpPlp/BsD7uogKIZYQPlbg1HZZzMTI1kcs6fvh2\nFln7G0uEsybI2QSnyfRZ5CKo9mNzgqEksFcMoP0BvC1AS77VoPFnqGqEmTCpVXz2ICSaj69gMy1D\nFZPZXE+w+2OTdWqJHyBiqkpdrwJ8e6j3/uMmmDk8oI3WAk0oW9/0C2BpiAQbKuobCv93GVTHzsvA\nuIYZDMtDwp8X7TbIxl1SS6WLo4RCJ8H5z8nglHSEvw6shGCB8QUZrSr9td5ZuWASbWY1NjjUk8cj\n2M622+UmUUusFj6eYVhjdlFmSg3whJV7S7naWSadhIFK2BBgOgKCYLd0RMASH78oh2Wr8nZg0FZb\n2GnsGXgYe84v6PozfWBg8U/nM1CXltzqUa4i5QY1YGQdMEFBqUciXpYg9vfPrjCOuJ2hw69w6HhK\njK+vERr4lIk2X0GpTxR9ljdyWDWcyIyqOkOJ3gJ+8eXFL8Sz8ne0/FZrdqUL9JTs9/7ptBWiCU/N\nVZ2hRHxIbRMXLE9c9TNBbtWWssJSh5I8VKX4TAvQRKo3bjU4FEhPYZd6GIs3LuVQ75aiYqx8BCqu\nlbVi/6cRF56QiDM96/ADLzVxysglTZxLqhZe0bf3T8zxmj+8dNckSUkbyHqFms0jmEesHDefs/yb\nTbvCC/YTkpBvAT8WWHxWeoGXzSG4lkjc4jjgvFLvsSu9YCeuLlhX7xyo+8ULrNwcxqcgnwR8j6sb\n2l3pBdxyBvJJEMyBD3EPit6vXqX+PJZ4kB4pdRh8LHhZeo6nvUyPf7PWGt4mmOj3Nj9guIcJ/gNu\n5vy5m96gbQAAAABJRU5ErkJggg==\n", | |
"prompt_number": 107, | |
"text": " ____________________________________\n \u2571 2 \n \u2571 ___ -erfinv (f) \n\u2572\u2571 \u2572\u2571 \u03c0 \u22c5f - 2\u22c5\u212f \u22c5erfinv(f) \n\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 ___ \n \u2572\u2571 \u03c0 " | |
} | |
], | |
"prompt_number": 107 | |
}, | |
{ | |
"cell_type": "code", | |
"collapsed": false, | |
"input": "# ratio for f = 0.9\nN(ratio.subs(f, 0.9))", | |
"language": "python", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"latex": "$$0.748808343784453$$", | |
"output_type": "pyout", | |
"png": "iVBORw0KGgoAAAANSUhEUgAAAKoAAAAPBAMAAACGiUnsAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEJmJdjLNVN0iZu+7\nq0QgoRR7AAAACXBIWXMAAA7EAAAOxAGVKw4bAAACsUlEQVQ4Ea2Uv2sTYRjHv3fXXJpLmhzWoeCQ\n2NJNajAOohUCFh0cLBXcxFQFF8E4dBGkxcWlxVBB0KWHFkVFbAcRsdbrJC42SEUohAYRFwetltb+\n8vy+z3PiP+BL8rm79/O9J3fvjwA79neDLTHWVSrxOADsPNITQ5zVt6+JVKm3amLUN0tv4rQY3lqg\nMDjvDzRFAxfQUWOvF0XRNo+jQAFtVYU4G/iIcSQ3aKlTBQwHmhYDtIYUBq+jFdXwpuBUANgsXQda\nZpAMkKoI1F0BbuFEFb+Yom5ZQ25E02KA4yGNwbmFmmpk6kgus9fhpx+4cxFOCGtEoK4X6MDXprXK\nFHV2FvmKpsUgdZ8FBXwq1cjVkV3jFdtnfou8bct3AoG6b3M4SCMjQM02HBBMq3HckE9kYKqyUecL\nyK7oVcjRrfG2pQ1WMVDXFt0LGHALsQbO8hIhoOahy1PB4qRMPPVEEYmfDPHnAs48WNWOBhWx27vq\nAyevMSI6caAQp2GMVXY5ZAaYw0QTov9VzYMDYKrePbXZFKjzui/NUjnPYg0c4gSYtBgHbshZIdgy\nIwR1vvh3BI5y7PpZ1avj6bRA3RN4m3xYvPRF8yw9DTANMUNwQwjYY68T1LkCkjpbo0A7WDXtI7Em\nUMcSZ6q3gfmq6FQNGd7ANIy5XIQbWgZoq8D5rZqr1zMrCxb3wNVGY3PBvNuiQJz1ndMSRD7mA9G5\nZWQ2JC3mQaOx9OqYQTldgb2sumUKdsVUTfwwxAzSNeCtQN0LDldtBrjBfqNDtK5rWgzfOKQgbO6w\nimpcx65+s8KTWnUL2UF4ZYE4fPExiUdIyvrbghdguKxpMUDOVCW4lx/XVKO97wMwxt3xnA690SJO\nd/VAIS7xnv8u2c49zVh/6nwXp8XAmd8OFOOHdwOiTan/3v4AQ1kYJE5RNsQAAAAASUVORK5CYII=\n", | |
"prompt_number": 108, | |
"text": "0.748808343784453" | |
} | |
], | |
"prompt_number": 108 | |
} | |
], | |
"metadata": {} | |
} | |
] | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment