Revision history [back]
 | 1 | initial version |
Your first (original) attempt can be fixed by making k and K symbolic, and using substitution:
import scipy.stats
n=20
p=0.05
binom_dist = scipy.stats.binom(n,p)
Bin=bar_chart([binom_dist.pmf(x) for x in range(n+1)],width=1)
T = RealDistribution('gaussian', 1)
#
var('k')
a = 1/(9*n-9*k)
b = 1/(9*k+9)
r = (k+1)*(1-p)/(n*p-k*p)
c = (1-b)*r^(1/3)
μ = 1 - a
σ = sqrt(b*r^(2/3) + a)
K = (c - μ)/σ
#
CP=plot(lambda k: T.distribution_function(K.subs(k=k)),(k,-1,21), rgbcolor=(0.8,0,0))
show(Bin+CP)
Note that in K.subs(k=k) the first k is the name of a keyword parameter (which will be interpreted by Sage as the name of a symbolic variable), and the second k is the argument which is passed into the lambda function.
Your second (original) attempt can be fixed by making k symbolic and doing the substitution inside distCPapprox:
import scipy.stats
n=20
p=0.05
binom_dist = scipy.stats.binom(n,p)
Bin=bar_chart([binom_dist.pmf(x) for x in range(n+1)],width=1)
T = RealDistribution('gaussian', 1)
#
var('k')
a = 1/(9*n-9*k)
b = 1/(9*k+9)
r = (k+1)*(1-p)/(n*p-k*p)
c = (1-b)*r^(1/3)
μ = 1 - a
σ = sqrt(b*r^(2/3) + a)
def distCPapprox(k): return T.distribution_function(((c - μ)/σ).subs(k=k))
#
CP=plot(distCPapprox,(k,-1,21), rgbcolor=(0.8,0,0))
show(Bin+CP)
Finally, the solution in your own comments is to create a function that takes a number and returns a number (avoiding symbolics entirely); personally I very much prefer that approach if symbolics are not needed.
