[SciPy-Dev] Stochastic calculus package
Maurizio Cipollina
maurizio.cipollina at gmail.com
Fri Jun 15 17:45:02 EDT 2018
Thanks Josef, point very well taken, the look is definitely old fashioned -
back from the good old times when a string wouldn’t be a String, and a list
wouldn’t be a List...
Maurizio
On Fri, 15 Jun 2018 at 20:39, <josef.pktd at gmail.com> wrote:
> Sounds good,
>
> One thing that immediately jumps into my eyes are the lower case class
> names.
> I think using pep-8 capital/camel case would make it look more modern.
> Only old fashioned packages like numpy and scipy still have the legacy
> spelling (when developers thought python has to look like matlab :)
>
> Josef
>
>
>
> On Fri, Jun 15, 2018 at 1:51 PM, Maurizio Cipollina
> <maurizio.cipollina at gmail.com> wrote:
> > Hi everybody,
> >
> > following up on this conversation of some time ago, a pre-release of the
> > package, renamed SdePy, is now available on pip
> > (https://pypi.org/project/sdepy ) and as a project on github
> > (https://github.com/sdepy/sdepy ), version 1.0.0rc3. A short quickstart
> > guide (https://sdepy.readthedocs.io/en/latest/intro.html#id2 ) is the
> > easiest path to see what it does and what it does not. Any feedback,
> > comments and suggestions are very welcome!
> >
> > Since last September, the package has evolved quite a bit and had its
> fair
> > share of refactoring, going through several 0.x unpublished versions
> before
> > stabilizing in its current form (marked as beta). Here are a few
> takeaways I
> > got from this thread, which I acted upon:
> >
> > - Integrating user defined SDEs turned out to be a main point of
> interest:
> > this ‘odeint’ like functionality has been considerably streamlined and
> > simplified - formerly, preset processes were mainly in focus, and
> getting a
> > new process up and running took a deep dive in a not-so-welcoming API,
> now
> > upgraded as well (the former ‘integrator’ has been superseded by a
> > ‘paths_generator’ parent class, and cooperating ‘integrator’ and ‘SDE’
> > classes).
> >
> > - The ‘kfunc’ beast was admittedly a hard sell: all undue dependencies
> of
> > the core modules on its machinery have been eliminated, it plays now the
> > role of an optional decorator for the shortcuts of frequently used stuff
> > (it’s fully functional though and, in my experience, quite handy).
> >
> > - The question “Do you REALLY need a ‘process’ class” clearly stood out
> in
> > your comments: after some experimentation, the main reason for answering
> in
> > the yes turned out, in my view, to be its seamless interoperability with
> the
> > rest of the package, hardly achievable without a dedicated design. Now
> > process instances can indeed be fed to integrators, both as
> time-dependent
> > and/or path-dependent SDE parameters, and as stochasticity sources: this
> > opened countless hacking opportunities when doing computations, and
> spoiling
> > them seemed a pity.
> >
> > In case any of you is willing to take a look, feel free to get in touch
> for
> > any problems/bugs/questions you might have
> >
> > Cheers
> > Maurizio
> >
> > On Mon, 11 Sep 2017 at 19:39, <josef.pktd at gmail.com> wrote:
> >>
> >>
> >>
> >> On Mon, Sep 11, 2017 at 12:35 PM, Maurizio Cipollina
> >> <maurizio.cipollina at gmail.com> wrote:
> >>>
> >>> Hi Ralf,
> >>>
> >>>
> >>> thanks for going through this, here are some follow up comments on your
> >>> points:
> >>>
> >>>
> >>>
> >>> On the two packages you found:
> >>>
> >>> Sdeint: provides several SDE integration algorithms but no specific
> >>> processes; does not cover processes with jumps; offers no explicit
> >>> provisions for handling solutions in a number of paths, which is what
> you
> >>> need since the integration output is random variables, not numbers.
> >>> SDE-higham is not a package or tool, but rather a collection of recipes
> >>> for a few textbook examples (no jumps either), in the form that an
> >>> interactive session might take (set globals, do calculations, plot and
> print
> >>> results).
> >>>
> >>>
> >>>
> >>> About your questions/concerns:
> >>>
> >>> ndarray subclassing: the integrator does not depend on the process
> class,
> >>> and process realizations may be easily tweaked to return an ndarray
> instead
> >>> of a process instance; however, a main point of the package is not
> just to
> >>> produce SDE solutions, but to allow the user to handle them (almost)
> >>> effortlessly once obtained: so one might possibly hide the process
> class
> >>> inside the testing suite, which depends on it, but hardly get rid of it
> >>> altogether. Indeed, I subclassed ndarray as safely as I could figure
> out,
> >>> but I will surely go through the scipy-dev archives and find out what
> >>> advised you against doing so.
> >>> odeint and scipy.integrate: my view is that odeint churns out numbers
> >>> dependant on numbers, while SDE integration churns out random variables
> >>> dependant on random variables (the stochasticity sources), so there is
> a
> >>> shift of paradigm and the two schemes do not fit smoothly one within
> the
> >>> other. The step by step integration of SDEs makes sense, and gets
> useful, if
> >>> you can easily control the inputs (eg. correlations among stochasticity
> >>> sources fed into different SDEs), generate the output in a decent
> number of
> >>> paths, and easily extract statistics of expressions dependent on it
> (hence
> >>> the flow ‘get a number of paths packaged as a process instance’ -> ‘do
> numpy
> >>> algebra with it’ -> ‘estimate probability distributions etc. using the
> >>> process class methods’ -> ‘scale it up with the montecarlo class in
> case the
> >>> needed paths do not fit into memory’).
> >>> Montecarlo simulations is a big topic indeed, that goes far beyond this
> >>> package; what is called the montecarlo class is in fact a tool to
> cumulate
> >>> statistical information extracted from a number of realizations of a
> random
> >>> variable, and is focused on making easy to interpret the output of SDE
> >>> integration.
> >>> Black-Scholes assists one more tests in a rather extensive testing
> suite,
> >>> and might be painlessly dropped.
> >>>
> >>>
> >>>
> >>> This said, I see your general feeling is that this whole subject is too
> >>> specialized to fit into scipy at this stage. I suggested otherwise, in
> my
> >>> proposal, because stochastic calculus is regarded as one of the
> pillars of
> >>> modern probability theory, and covering it into scipy might be
> relevant to
> >>> lots of mathematicians and physicists, both practitioners and
> students, that
> >>> might otherwise hesitate to adopt, and build dependencies upon, a
> standalone
> >>> package. In case you reconsider let me know, I’d be happy to help (now
> or in
> >>> the future…).
> >>>
> >>>
> >>>
> >>> If you have further questions, let me know as well.
> >>>
> >>> Cheers
> >>>
> >>> Maurizio
> >>
> >>
> >>
> >> To partially follow up, similar to Ralph's view.
> >>
> >> I think it would be good to have this as a separate package, or at least
> >> the code on github.
> >> With distributions like conda it is now much easier to install
> additional
> >> packages, once they are established enough and are included in
> distributions
> >> or easily pip installable.
> >>
> >> The main reason that I think scipy is currently not the appropriate
> >> package for it is that the design and review would require a lot of
> work. I
> >> guess that there are not many maintainers and reviewers that are
> familiar
> >> with this area. I also know only of applications in finance and it is
> not
> >> clear whether or which parts would be of more general interest.
> >>
> >> In my opinion, some core tools for SDE and jump diffusions would fit in
> >> scipy. But the application support will not or might not fit well in a
> >> general purpose scientific numerical library. This would be similar to
> other
> >> areas where scipy provides the core computational tools, but
> applications
> >> are supported by other packages.
> >>
> >> For example, the process class might be similar in magnitude to scipy
> >> stats distributions or signal statespace classes, but I have no idea
> what
> >> design for it would fit in scipy.
> >> Some of the time aggregation properties sound similar to the
> corresponding
> >> pandas functions, and it is not clear whether users would want to use a
> >> separate class or just use pandas instead.
> >>
> >> Similarly, the kfunc class sounds like something that doesn't have a
> >> similar pattern in scipy.
> >>
> >> Design changes can be made more easily in a standalone package, while
> all
> >> refactorings of classes in scipy are difficult issues because of the
> >> backwards compatibility policy which requires that the design should be
> >> largely settled before a merge.
> >>
> >> Josef
> >>
> >> https://groups.google.com/d/msg/pystatsmodels/xsttiEiyqAg/maR6n4EeAgAJ
> >> (stackoverflow question has been deleted "Is there a library out there
> to
> >> calibrate an Ornstein-Uhlenbeck model?")
> >>
> >>
> >>
> https://github.com/statsmodels/statsmodels/blob/master/statsmodels/sandbox/tsa/diffusion.py
> >> and diffusion2.py
> >>
> >>
> >>
> >>>
> >>>
> >>>
> >>>
> >>> On 9 September 2017 at 03:31, Ralf Gommers <ralf.gommers at gmail.com>
> >>> wrote:
> >>>>
> >>>> Hi Maurizio,
> >>>>
> >>>>
> >>>> On Fri, Sep 8, 2017 at 7:55 PM, Maurizio Cipollina
> >>>> <maurizio.cipollina at gmail.com> wrote:
> >>>>>
> >>>>> Hi everybody, and first of all thanks for the great job you have been
> >>>>> doing, on which I relied *A LOT* over time.
> >>>>>
> >>>>> Stochastic calculus and Stochastic Differential Equations are a
> >>>>> significant branch of mathematics that to my understanding is
> >>>>> underrepresented in scipy, as well as in the open source python
> stack at
> >>>>> large. I may well be missing something, but as I can see one can
> find out
> >>>>> there many recipies or collections of recipies, and a number of
> >>>>> packages/projects skewed towards finance, but no structured
> framework for
> >>>>> generating and handling stochastic processes and numerically solving
> SDEs.
> >>>>
> >>>>
> >>>> I've never used these, but a quick search turns up two packages:
> >>>> https://pypi.python.org/pypi/sdeint
> >>>> https://github.com/alu042/SDE-higham
> >>>> Not very mature probably, but would be good to compare your code with
> >>>> those.
> >>>>
> >>>>>
> >>>>> This is a pity, since scipy/numpy provide efficen implementations of
> >>>>> all the basic ingredients needed to do the trick. Out of
> professional need,
> >>>>> and prompted by this gap, I have developed a general purpose package
> of
> >>>>> which I own the copyright, and which I would be happy to release
> under the
> >>>>> BSD license as part of scipy, in case there is a consensus among the
> scipy
> >>>>> community that it makes sense to do so (the package has no
> dependencies
> >>>>> beyond standard library, numpy and scipy, and might probably fit
> best as a
> >>>>> subpackage of scipy.stats).
> >>>>>
> >>>>> Before setting up a PR, I would be happy to hear your thoughts: I
> have
> >>>>> pasted below the package docstring, that should give an overview of
> its
> >>>>> scope and limitations. Please note that some of the realized
> processes stray
> >>>>> into mathematical finance territory, somehow inevitably since
> finance is one
> >>>>> major field of application of stochastic calculus, but the focus is
> the
> >>>>> latter, not the former, at least in my intentions.
> >>>>
> >>>>
> >>>> First thought: this looks like a useful addition to the scientific
> >>>> Python ecosystem, but is probably too specialized for SciPy (at least
> at
> >>>> this stage).
> >>>>
> >>>> Based on your docstring below, here are some questions/concerns:
> >>>> 1. We definitely don't want a new ndarray subclass (there's a pretty
> >>>> strong consensus at this point that subclassing ndarray is usually
> not a
> >>>> good solution), so your `process` doesn't sound attractive.
> >>>> 2. The `integrator` class sounds a lot like odeint, so maybe this
> would
> >>>> fit better with scipy.integrate?
> >>>> 3. Monte-Carlo simulation is a big and fairly technical topic. There
> are
> >>>> whole packages (e.g. PyMC3, emcee) dedicated to this. It may not fit
> well
> >>>> with SciPy.
> >>>> 4. Things specific to finance like Black-Scholes and put/call options
> >>>> are not a good fit.
> >>>> 5. Maintainers for this added code. The ODE integrators in
> >>>> scipy.integrate already suffer from lack of a dedicated maintainer,
> SDEs are
> >>>> more specialized so are likely to not be very well-maintained by
> other devs
> >>>> than you.
> >>>>
> >>>> Without seeing your code it's hard to say for sure, but it looks to me
> >>>> like it would be better to release your code as a standalone package.
> >>>>
> >>>> Cheers,
> >>>> Ralf
> >>>>
> >>>>
> >>>>>
> >>>>> Thanks in advance for taking your time to go through this, and for
> your
> >>>>> comments and suggestions.
> >>>>> Maurizio
> >>>>>
> >>>>> “””
> >>>>> ===========================================================
> >>>>> Stochastic - Monte Carlo simulation of stochastic processes
> >>>>> ===========================================================
> >>>>>
> >>>>> This package provides:
> >>>>>
> >>>>> 1. A `process` class, a subclass of `numpy.ndarray` representing
> >>>>> a sequence of values in time, realized in one or several paths.
> >>>>> Algebraic manipulations and ufunc computations are supported for
> >>>>> instances that share the same timeline and comply with numpy
> >>>>> broadcasting rules. Interpolation along the timeline is supported
> >>>>> via callability of `process` instances. Process-specific
> >>>>> functionalities,
> >>>>> such as averaging and indexing along time or across paths, are
> >>>>> delegated
> >>>>> to process-specific methods, attributes and properties (no
> >>>>> overriding
> >>>>> of `numpy.ndarray` operations).
> >>>>>
> >>>>> 2. The `source` class and its subclasses, stochasticity sources
> >>>>> providing
> >>>>> numerical realizations of the differentials commonly found
> >>>>> in stochastic differential equations (SDEs), with or without
> >>>>> memory of formerly invoked realizations.
> >>>>>
> >>>>> 3. The `integrator` class, a general framework for numerical SDE
> >>>>> integration, intended for subclassing, to assist the definition
> of
> >>>>> step by
> >>>>> step integration algorithms and the computation of numerical
> >>>>> realizations of
> >>>>> stochastic processes.
> >>>>>
> >>>>> 4. Several objects providing realizations of specific stochastic
> >>>>> processes
> >>>>> along a given timeline and across a requested number of paths.
> >>>>> Realizations are obtained via explicit formulae, in case the
> >>>>> relevant SDE
> >>>>> has a closed form solution, or otherwise as `integrator`
> subclasses
> >>>>> performing
> >>>>> Euler-Maruyama numerical integration.
> >>>>>
> >>>>> 5. A `montecarlo` class, as an aid to cumulate the results of
> several
> >>>>> Monte Carlo simulations of a given stochastic variable, and to
> >>>>> extract
> >>>>> summary estimates for its probability distribution function and
> >>>>> statistics, thus supporting simulations across a number of paths
> >>>>> that exceeds
> >>>>> the maximum allowed by available memory.
> >>>>>
> >>>>> 6. Several analytical results relating to the supported stochastic
> >>>>> processes, as a general reference and for testing purpouses.
> >>>>>
> >>>>> For all sources and realizations, process values can take any shape,
> >>>>> scalar or multidimensional. Correlated multivariate stochasticity
> >>>>> sources are
> >>>>> supported. Time-varying process parameters (correlations, intensity
> of
> >>>>> Poisson
> >>>>> processes, volatilities etc.) are allowed whenever applicable.
> >>>>> `process` instances act as valid stochasticity source realizations
> (as
> >>>>> does
> >>>>> any callable object complying with a `source` protocol), and may be
> >>>>> passed as a source specification when computing the realization of a
> >>>>> given
> >>>>> process.
> >>>>>
> >>>>> Computations are fully vectorized across paths, providing an
> efficient
> >>>>> infrastructure for simulating a large number of process realizations.
> >>>>> Less so, for large number of time steps: integrating 1000 time steps
> >>>>> across 100000 paths takes seconds, one million time steps across
> >>>>> 100 paths takes minutes.
> >>>>>
> >>>>> General infrastructure
> >>>>> ======================
> >>>>> .. autosummary::
> >>>>> :toctree: generated/
> >>>>>
> >>>>> process -- Array representation of a process (a subclass of
> >>>>> -- numpy.ndarray).
> >>>>> kfunc -- Base class for functions with managed keyword
> >>>>> arguments.
> >>>>> source -- Base class for stochasticity sources.
> >>>>> true_source -- Base class for stochasticity sources with memory.
> >>>>> integrator -- General framework for numerical SDE integration.
> >>>>> montecarlo -- Summmary statistics of Monte Carlo simulations.
> >>>>>
> >>>>>
> >>>>> Stochasticity sources
> >>>>> =====================
> >>>>> .. autosummary::
> >>>>> :toctree: generated/
> >>>>>
> >>>>> wiener -- dW, a source of standard Wiener process
> >>>>> (Brownian
> >>>>> -- motion) increments.
> >>>>> symmetric_wiener -- as above, with symmetric paths (averages
> >>>>> exactly 0).
> >>>>> true_wiener -- dW, a source of standard Wiener process
> >>>>> (Brownian
> >>>>> -- motion) increments, with memory.
> >>>>> poisson -- dN, a source of Poisson process increments.
> >>>>> compound_poisson -- dJ, a source of compound Poisson process
> >>>>> increments.
> >>>>> -- (jumps distributed in time as a Poisson
> >>>>> process,
> >>>>> -- and in size as a given `scipy.stats`
> >>>>> distribution).
> >>>>> compound_poisson_rv -- Preset jump size distributions for compound
> >>>>> Poisson
> >>>>> -- process increments.
> >>>>>
> >>>>>
> >>>>> Stochastic process realizations
> >>>>> ===============================
> >>>>> .. autosummary::
> >>>>> :toctree: generated/
> >>>>>
> >>>>> const_wiener_process -- Wiener process (Brownian motion),
> >>>>> with
> >>>>> -- time-independent parameters.
> >>>>> const_lognorm_process -- Lognormal process, with
> >>>>> time-independent
> >>>>> -- parameters.
> >>>>> wiener_process -- Wiener process (Brownian motion).
> >>>>> lognorm_process -- Lognormal process.
> >>>>> ornstein_uhlenbeck_process -- Ornstein-Uhlenbeck process
> >>>>> (mean-reverting
> >>>>> -- Brownian motion).
> >>>>> hull_white_process -- F-factors Hull-White process (sum
> >>>>> of F
> >>>>> -- correlated mean-reverting
> >>>>> Brownian
> >>>>> -- motions).
> >>>>> hull_white_1factor_process -- 1-factor Hull-White process (F=1
> >>>>> Hull-White
> >>>>> -- process with F-index collapsed
> >>>>> to a
> >>>>> -- scalar).
> >>>>> cox_ingersoll_ross_process -- Cox-Ingersoll-Ross mean-reverting
> >>>>> process.
> >>>>> full_heston_process -- Heston stochastic volatility
> >>>>> process
> >>>>> -- (returns both process and
> >>>>> volatility).
> >>>>> heston_process -- Heston stochastic volatility
> >>>>> process
> >>>>> -- (stores and returns process
> >>>>> only).
> >>>>> jump_diffusion_process -- Jump-diffusion process (lognormal
> >>>>> process
> >>>>> -- with compound Poisson jumps).
> >>>>> merton_jump_diffusion_process -- Merton jump-diffusion process
> >>>>> -- (jump-diffusion process with
> >>>>> normal
> >>>>> -- jump size distribution).
> >>>>> kou_jump_diffusion_process -- Kou jump-diffusion process
> >>>>> (jump-diffusion
> >>>>> -- process with double
> exponential
> >>>>> -- jump size distribution).
> >>>>>
> >>>>> Shortcuts::
> >>>>>
> >>>>> lognorm -- lognorm_process
> >>>>> oruh -- ornstein_uhlenbeck_process
> >>>>> hwp -- hull_white_process
> >>>>> hw1f -- hull_white_1factor_process
> >>>>> cir -- cox_ingersoll_ross_process
> >>>>> heston -- heston_process
> >>>>> mjd -- merton_jump_diffusion_process
> >>>>> kou -- kou_jump_diffusion_process
> >>>>>
> >>>>>
> >>>>> Analytical results
> >>>>> ==================
> >>>>> Exact statistics for the realized stochastic processes
> >>>>> are listed below, limited to the case of constant process parameters
> >>>>> and with
> >>>>> some further limitations to the parameters' domains.
> >>>>> Function arguments are consistent with those of the corresponding
> >>>>> processes.
> >>>>> Suffixes `pdf`, `cdf` and `chf` stand respectively for probability
> >>>>> distribution
> >>>>> function, cumulative probability distribution function, and
> >>>>> characteristic
> >>>>> function. Black-Scholes formulae for the valuation of call and put
> >>>>> options have been
> >>>>> included (prefixed with `bs` below).
> >>>>>
> >>>>> .. autosummary::
> >>>>> :toctree: generated/
> >>>>>
> >>>>> wiener_mean
> >>>>> wiener_var
> >>>>> wiener_std
> >>>>> wiener_pdf
> >>>>> wiener_cdf
> >>>>> wiener_chf
> >>>>>
> >>>>> lognorm_mean
> >>>>> lognorm_var
> >>>>> lognorm_std
> >>>>> lognorm_pdf
> >>>>> lognorm_cdf
> >>>>> lognorm_log_chf
> >>>>>
> >>>>> oruh_mean
> >>>>> oruh_var
> >>>>> oruh_std
> >>>>> oruh_pdf
> >>>>> oruh_cdf
> >>>>>
> >>>>> hw2f_mean
> >>>>> hw2f_var
> >>>>> hw2f_std
> >>>>> hw2f_pdf
> >>>>> hw2f_cdf
> >>>>>
> >>>>> cir_mean
> >>>>> cir_var
> >>>>> cir_std
> >>>>> cir_pdf
> >>>>>
> >>>>> heston_log_mean
> >>>>> heston_log_var
> >>>>> heston_log_std
> >>>>> heston_log_pdf
> >>>>> heston_log_chf
> >>>>>
> >>>>> mjd_log_pdf
> >>>>> mjd_log_chf
> >>>>>
> >>>>> kou_mean
> >>>>> kou_log_pdf
> >>>>> kou_log_chf
> >>>>>
> >>>>> bsd1d2
> >>>>> bscall
> >>>>> bscall_delta
> >>>>> bsput
> >>>>> bsput_delta
> >>>>>
> >>>>> Notes
> >>>>> =====
> >>>>> To the benefit of interactive and notebook sessions, and as an aid to
> >>>>> fluently
> >>>>> freeze or vary the plurality of parameters that define each
> stochastic
> >>>>> process,
> >>>>> all sources, process realizations and analytical results are objects
> >>>>> with managed keywords (subclasses of `kfunc`): if the corresponding
> >>>>> classes
> >>>>> are instantiated with both variables (non keyword arguments) and
> >>>>> parameters
> >>>>> (keyword arguments), they behave as functions returning the
> >>>>> computations'
> >>>>> result; if called with parameters only, return an instance that
> stores
> >>>>> the set
> >>>>> parameters, and exposes the same calling pattern (call with variables
> >>>>> and optionally with parameters to get results, call with parameters
> >>>>> only
> >>>>> to get a new instance storing modified parameters).
> >>>>> Parameters that are not specified fall back to the class defaults, if
> >>>>> calling
> >>>>> the class, or to the instance's stored values, if calling an
> instance.
> >>>>> """
> >>>>>
> >>>>> _______________________________________________
> >>>>> SciPy-Dev mailing list
> >>>>> SciPy-Dev at python.org
> >>>>> https://mail.python.org/mailman/listinfo/scipy-dev
> >>>>>
> >>>>
> >>>>
> >>>> _______________________________________________
> >>>> SciPy-Dev mailing list
> >>>> SciPy-Dev at python.org
> >>>> https://mail.python.org/mailman/listinfo/scipy-dev
> >>>>
> >>>
> >>>
> >>> _______________________________________________
> >>> SciPy-Dev mailing list
> >>> SciPy-Dev at python.org
> >>> https://mail.python.org/mailman/listinfo/scipy-dev
> >>>
> >>
> >> _______________________________________________
> >> SciPy-Dev mailing list
> >> SciPy-Dev at python.org
> >> https://mail.python.org/mailman/listinfo/scipy-dev
> >
> >
> > _______________________________________________
> > SciPy-Dev mailing list
> > SciPy-Dev at python.org
> > https://mail.python.org/mailman/listinfo/scipy-dev
> >
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at python.org
> https://mail.python.org/mailman/listinfo/scipy-dev
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20180615/5c02026f/attachment-0001.html>
More information about the SciPy-Dev
mailing list