DEP: Depcrecate LTI classes in signal module

I would like to propose the deprecation of the “linear time invariant system base class” scipy.signal._ltisys.LinearTimeInvariant and its children. This post was motivated by the current discussions about deprecations for SciPy 2.0. Issue #6137 (opened in 2016) provides some further insights on this topic.

The object-oriented design of the classes in file scipy/signal/_ltisys.py relies
heavily on inheritance, i.e.:

class LinearTimeInvariant:
    """Linear time invariant system base class. """
    
class lti(LinearTimeInvariant):
    """Continuous-time linear time invariant system base class. """
    
class dlti(LinearTimeInvariant):
    """Discrete-time linear time invariant system base class. """
    
class TransferFunction(LinearTimeInvariant):
    """Linear Time Invariant system class in transfer function form. """
    
class TransferFunctionContinuous(TransferFunction, lti):
    """Continuous-time Linear Time Invariant system in transfer function form. """
    
class TransferFunctionDiscrete(TransferFunction, dlti):
    """Discrete-time Linear Time Invariant system in transfer function form. """
    
class ZerosPolesGain(LinearTimeInvariant):
    """Linear Time Invariant system class in zeros, poles, gain form. """
    
class ZerosPolesGainContinuous(ZerosPolesGain, lti):
    """Continuous-time Linear Time Invariant system in zeros, poles, gain form. """
    
class ZerosPolesGainDiscrete(ZerosPolesGain, dlti):
    """Discrete-time Linear Time Invariant system in zeros, poles, gain form. """
    
class StateSpace(LinearTimeInvariant):
    """Linear Time Invariant system in state-space form. """
    
class StateSpaceContinuous(StateSpace, lti):
    """Continuous-time Linear Time Invariant system in state-space form """
    
class StateSpaceDiscrete(StateSpace, dlti):
    """Discrete-time Linear Time Invariant system in state-space form """

The inheritance reduces some code duplication at the price of adding significant code
complexity and maintenance costs. E.g.:

  1. The classes LinearTimeInvariant, lti, dlti, TransferFunction, TransferFunctionContinuous, TransferFunctionDiscrete, ZerosPolesGain, ZerosPolesGainContinuous, and ZerosPolesGainDiscrete do not contain algorithms. Their methods just wrap other public functions.
  2. Instantiating a StateSpace produces either an instance of StateSpaceContinuous or StateSpaceDiscrete depending on the parameters. The same holds for ZerosPolesGain. This does not mix too well with the Python typing system and the new Array API backends.
  3. The functionality is inconsistent. E.g., StateSpace implements a multiplication operator (concatenation), but TransferFunction and ZerosPolesGain do not.
  4. The functions lsim, dlsim, impulse, dimpulse, dfreqresp, freqresp, step, and dstep accept instances from the classes above, but do not automatically distinguish between continuous-time and discrete-time systems.

Progress of new features and bug fixing of these classes is quite slow. In my opinion, a key reason for the few new contributions is the relatively steep learning curve for understanding the code.

Further reasons for not keeping these classes:

  • In theory these classes allow easily switching system representations. In practice, numerics frequently foil this. IMO, it is better to make the user call explicit conversion functions, which ideally would document numerical pitfalls.
  • It would not be too hard to implement system composability functionality by adding +, -, *, /, and @ operators, where applicable. This would make these classes much more versatile. But again, the numerics of almost stable systems are a rich source for potential issues.
  • Practical reasoning about continuous-time systems is quite different from reasoning about discrete-time systems. Hence, I think that mixing discrete-time and continuous-time systems in the same class is a little bit confusing from a user’s perspective: A newcomer might only be well-versed in one of them.

For that reason I would propose deprecating the classes listed above. Furthermore,
to depecreate that the functions lsim, dlsim, impulse, dimpulse, step, dstep, dfreqresp, freqresp, step, dstep, ss2tf, ss2zpk, cont2discrete, and place_poles accept instances of those classes as an argument.

Outlook

Having a modernized set of LTI classes is a worthwhile goal for SciPy, IMO.
If such an ambitious project is pursued, the following features would be desirable:

  • There would have to be well-documented mechanisms to deal with the most common numerical problems. E.g., detecting and eliminating (almost) equal zeros and poles in transfer functions or detect (almost) unstable state space systems.
  • The functionality of the functions lsim, dlsim, impulse, dimpulse, dfreqresp, freqresp, step, dstep, ss2tf, ss2zpk, cont2discrete, place_poles should become class methods.
  • System composability by adding +, -, *, /, and @ operators, where applicable. Also, adding functionality for designing (linear) controllers like PID or LQG would be a boon.
  • It would be great to have a pluggable arbitrary-precision floating-point arithmetic backend like mpmath. E.g., this would allow calculating eigenvalues with arbitrary
    accuracy.

I would love to read feedback!

What (roughly) would the part of the release notes look like which directs users of the deprecated functionality to alternatives?

The release notes could look something like this:

The user-facing classes StateSpace, TransferFunction, ZerosPolesGain, lti, dlti are deprecated. Their children, i.e., LinearTimeInvariant, TransferFunctionContinuous, TransferFunctionDiscrete, ZerosPolesGainContinuous,
ZerosPolesGainDiscrete, StateSpaceContinuous, and StateSpaceDiscrete as well as their abstract base class LinearTimeInvariant are deprecated as well.

All their functionality, (i.e., their methods) are already present as separate functions. Use the LTI representations functions to convert between the system representations.

Instead of passing instances of those classes to the functions lsim, dlsim, impulse, dimpulse, dfreqresp, freqresp, step, dstep, ss2tf, ss2zpk, cont2discrete, place_poles, do (this already works):

  • StateSpace → the tuple (A, B, C, D) made up of made four 2d arrays
  • ZerosPolesGain → the tuple (zeros, poles, gain) made up of two 1d arrays and a scalar
  • TransferFunction → the tuple (numerator, denominator) made up of two 1d arrays

FWIW; The current class hierarchy is pretty straightforward to write accurate stubs for in scipy-stubs (src). It would indeed be simpler if there’d just be one way to do things, so I’m +1 on removing these. But I’d also be +1 if the proposal was to depracate the base tuple-based LTI representations in favor of the (more explicit) class-based interface.

But implementation complexity aside, I think that the ergonomics from the user’s perspective should take priority here. That said, I personally don’t use this LTI functionality in any of my projects, so I’m probably not the right person to ask.

I personally do not really use them much either. I agree with you that a more explicit class-based interface is preferable. IMO, the current class hierarchy is kind of upside-down, since it relies (more or less) on implicit system representation conversions, which can mask numerical problems.

I did reflect a bit on how to design a more intuitive interface. Currently, what seems most promising to me is implementing each class independently of the others, with explicit methods for converting to different system representations. E.g., this would make it straightforward to implement different freqresp methods for different representations accounting for numerics..

What made my thinking process stall:

  • I wanted to read a bit about fixed-point filter design to see if it is possible to incorporate some of those ideas. I have not really found the time yet, nor a good book on the topic.
  • What the best way is to help users detect potential numerical problems.
  • The scope of future functionality (e.g., controller design, etc.).

TL;DR in my opinion the options are: either we sweep a large portion of signal away or we commit to higher quality and do it properly. Anything in between is user discomfort for very little benefit in my opinion. Currently, no serious work can be done with the existing tooling.

Quality wise I agree that LTI and co. are not of contemporary quality. They implement the “naive” algorithms and not production grade versions, in the reasonable meaning of “production” hence I would support the idea of the clean-up.

However things are a bit more complicated than that when it comes to filtering since the rest of the signal module operates through these models as currency. Since matlab times, signal and control toolboxes always have an annoying compatibility issues, say FIR coefficients are used in reversed order and so on. But in general the objects are recycled, that is to say, they all use TF, SS, ZPK models.

We can argue that these implementations can be converted into hidden private API but that is also quite awkward to have. In the python-control and harold packages you can see that the edge cases are quite many and in fact they are the most wasted-effort parts of these packages. There are lots of edge cases that needs attention, say multiple by scalar, MIMO details, arithmetic with numpy arrays of size-matching systems, say 2 x 4 matrix multiplying 4 input 3 output LTI object and so on. So either way we need proper (no pun intended) LTI models for the rest.

Indeed fixing these functions are quite nontrivial. Even deciding on the time span of a step response is jumping through loads of hoops. As an example, this is something folks moved from harold to control python-control/control/timeresp.py at 7d5019e0b31f7f9f1099808d15cedfbaf01a0a56 · python-control/python-control · GitHub which is not something you arrive at by simple manipulations of system objects. Same with bode and nyquist plots. There are much better ways to compute the frequency response over the “naive” C @ solve(iw-A, B) + D calls which are numerically not so robust and very slow per w.

You typically don’t need and in fact refrain from doing so. I did it in harold but it is not solving any problems but in fact creating more as sometimes numerics disturb the results leading to user surprise though it is probably doing the right thing. This is best left to the user to run minimal realizations over the LTI models.

It does not have to be class methods but can stay functions. Their quality needs to be improved. OOP does not solve the underlying issue. And on top of that you have to implement three copies of per model. Standalone functions can accept all three and act accordingly.

Composability is similar to above, not very easy to do in a straightforward fashion especially ill-posed block diagrams cause a lot of issues, say, nonzero feedthrough cases and so on. control package might be a better target for control synthesis though we have a pole placement in SciPy already. Again a sign of inconsistent API detailing and legacy.

This does not matter too much since as the array size increases mpmath also suffers from the same issues as eigenvalues for n>= 5 is inaccessible in closed form. The pole locations are typically good enough for virtually all use cases, ala “if it is small enough it is a zero eigenvalue”. The issue is mostly conversion precision from one model to the next and discretization care.

Back in the day, I was thinking about proposing SLICOT fortran library (now I have the C-translated version) but it is a bit too specialized to control purposes (it also has lots of unnecessary linalg and sysid parts). We can scavenge harold for parts as I don’t have anything planned in the near future for it due to time constraints. Or leave all functionality to control and work on expectation management for the users.

I think the class-based implementations are causing most of the issues, as the saying goes, “composition over inheritance” is not followed and everything is leaking internal abstractions which is very common for the early 2000s code. We learned a lot since then how to handle classes and these models should stand independent without insisting on all models have all methods from the base class etc. They have too many edge cases to have common representations.

Thank you Ilhan for putting the focus on the algorithmic side—the term “cleanup” hits the nail on pretty much on the head.

I fully agree. As you know, the things that age much faster than software are decisions on software design.

I trust your judgment on that. Though, in my perception, a not insignificant number of issues stems from users not being aware that they are dealing with numerically ill-posed systems. Some advice in the documentation on how to recognize those would go a long way.

Regarding the arbitrary-precision floating-point arithmetic, the main use case I had in mind was filter design: Using mpmath for getting a ground-truth when simulating designs would have been quite useful to me in the past. I do see though that it would be quite ambitious to integrate.

I was not aware harold existed. Though it seems to have less functionality than python-control, it user interface seems more accessible to me.

It is a good point to ask which functionality should be improved and which should be left to python-control. I.e., it would make sense to take a little time to a high-level roadmap for the LTI systems before moving on. A starting point for this could be:

  1. Deprecate the LTI classes (as proposed above) and point the users to python-control / harold (as noted, the classes do not contain any algorithms)
  2. Remove the internal usage of those classes (it’s mostly mechanical work).
  3. Clarify the documentation of freqresp, dfreqresp that they are only a thin wrappers of the functions freqs, freqs_zpk as well as to freqz and freqz_zpk. Furthermore, clarify for the other functions in the LTI systems section, which internal system representation they use.
  4. The hard part: Discuss if there are other LTI systems functions, which should be marked as deprecated or legacy (I do not have strong opinions on this topic).
  5. The harder part: Improve the quality of the remaining functions and decide what functionality is missing.