I don't understand how to use cumulative distribution functions in sagemath. For example, I'd like to define a function N(x) as N(x)=1-F(x) where F is a cumulative log-normal distribution However, if I try
W = RealDistribution('lognormal',[1.5,.6])
N(x) = W.cum_distribution_function(x)
it does not seem to work at all
W.cum_distribution_function
is not a symbolic function. What you can do is creating N as an alias:
N = W.cum_distribution_function
print N(4)
plot(N)
But what if I need to use it as/in a symbolic function?
Using the actual function I get a fine plot via:
sage: CDFlogn(x,nu,sigma) = 1/2*(1+erf((log(x)-nu)/(sigma*sqrt(2))))
sage: plot(CDFlogn(x,1.5,.6),0,10)
Launched png viewer for Graphics object consisting of 1 graphics primitive
sage: CDFlogn(x,1.5,.6)
1/2*erf(-0.833333333333333*sqrt(2)*(-log(x) + 1.50000000000000)) + 1/2
sage: CDFlogn(4,1.5,.6)
1/2*erf(-0.833333333333333*sqrt(2)*(-log(4) + 1.50000000000000)) + 1/2
sage: CDFlogn(4,1.5,.6).n()
0.424846794957399
For example: https://en.wikipedia.org/wiki/Log-nor...
Use a Python function, either using def, or using lambda.
This will let you evaluate the function, and plot it.
Here, using lambda:
sage: W = RealDistribution('lognormal',[1.5, .6])
sage: def N(x): return W.cum_distribution_function(x)
sage: N(4)
0.4248467949573992
sage: plot(N, (x, 0, 10))
Here, using def:
sage: W = RealDistribution('lognormal',[1.5, .6])
sage: N = lambda x: W.cum_distribution_function(x)
sage: N(4)
0.4248467949573992
sage: plot(N, (x, 0, 10))
The down side is you need to do that for every variation based on N.
For instance, to plot 1 - N,
sage: NN = lambda x: 1 - N(x)
sage: plot (NN, (x, 0, 10))
Of course, to eg differentiate it, you would really need a symbolic function.
Maybe someone knows how to turn it into a symbolic function.
Maybe it's worth opening a Sage trac ticket to make cumulative distribution functions symbolic.
Thanks. However, this does not seem to allow me to use N(x) in an expression as I would do with any other (symbolic) function. For example, I'm not even able to 'plot(1-N,(x,0,10))' This is unfortunate, as I wanted to translate into sagemath a simulation I did in R, where I can use "proper" functions, such as 'plnorm' and 'dlnorm'
Here is a shortened version of what I did in R, and would like to replicate in sagemath:
x=(0:80)/10
meanl=5/3
sdl=2/3
n <- function(z) { 1-plnorm(z,meanlog=meanl,sdlog=sdl) }
b <- function(z) dlnorm(z, meanlog=meanl, sdlog=sdl)
BB <- function(z) integrate(function(w) w * b(w), lower = z, upper = Inf)$value
B <- Vectorize(BB, vectorize.args='z')
h <- function(x) 40*exp(-x)
plot(x,h(x)+x,ylim=c(0,10),type="l")
lines(x,5-x-n(h(x))*h(x)+B(h(x)),xlim=c(0,max(x)),col=3)
Asked: 2016-11-03 23:57:25 +0100
Seen: 1,663 times
Last updated: Nov 04 '16
Copyright Sage, 2010. Some rights reserved under creative commons license. Content on this site is licensed under a Creative Commons Attribution Share Alike 3.0 license.