by Richard Saumarez
You may wonder why a medic is writing a post on control theory in climate.
As Professor Curry asked me to give some biographical detail, I should explain that after medical school, I did a PhD in biomedical engineering, which before BME became an academic heavy industry, was in an electrical engineering department. My clinical research interest developed around the study of sudden cardiac death – people with genetic heart diseases who drop dead from a fatal rhythm disturbance of the heart. I soon realised that a purely clinical approach had to be augmented by quantitative understanding how basic cellular and physiological abnormalities in disease create the conditions for a fatal arrhythmia to develop; the problem being that measurements in an animal model may be very difficult to reproduce safely in humans. In clinical cardiac electrophysiology, we routinely stimulate the heart and make recordings from within it (it is surprisingly safe) and signal processing applied to intra-cardiac signals showed there was much more information in them than had been previously recognised. This allowed deductions about the mechanisms that precipitate ventricular fibrillation and prediction of sudden death in some patient groups. The link between basic physiology and man was established through mathematical models of the relevant cardiac electrophysiology so one can investigate the influence of basic pathological effects, reducing intercellular connectivity, changing ionic currents in cardiac cell membranes etc. on the measurements that one can make in man.
The use of engineering techniques in a clinical problem –sudden death – has taught me some of the pitfalls (by falling into them regularly) in interpreting data from one system in another and that carelessness in representing the basic physical nature of the system could lead to major conceptual problems. In physiology and medicine, the devil is the detail and a mathematical approach that fails to recognise these details may produce absurd results. Following the debate on climate feedback, I realised that there were similar problems: the physical concepts of climate and a formal approach to control and systems theory didn’t quite mesh.
Remembering AP Herbert’s famous dictum: “If nobody were to open their mouths without knowing exactly what they were talking about, a Deadly Hush would fall upon the World”, I have written a review of signal/control theory that is relevant to dynamic climate feedback. It is in two sections. This is a verbal, qualitative description of the problem, and my conclusions, of the argument that is developed at a simple mathematical level in the attached PDF [climate_feedback].
Physical assumptions in SB2011
I start with the equation used by Spencer & Braswell and Dessler to describe feedback. This is a first order differential equation. When you write a differential equation to describe something, in this case feedback, you have made a definite statement about how you think the systems works physically. In this case, the equation describes a simple physical system that contains an energy storage element and an energy dissipation element. Examples of this are: a resistor capacitor network, a spring/damper system reacting to a force or the concentration of injected dye passing through a mixing vessel. All these systems are described by the same equation and are said to be analogous.
|
I then discuss the process of “Convolution” and “deconvolution”. We are trying to identify the behaviour of a system, short-term climate feedback, by looking at goes in (net Flux) and what comes out (Temperature). A system can be defined in a number of ways, one of which is the impulse response (IR). This is the output of the system when its input is an impulse (an abstraction of an infinite magnitude, infinitely short pulse). If one knows the IR, because one can calculate it theoretically, it is straightforward to calculate the output from the IR and any input signal. This is called convolution and represents a laborious process of treating the input signal as a set of impulses, at infinitely small time steps, weighted by the amplitude at that instant, applying the impulse to each one and integrating them. Although convolution is a central analytical concept in systems and control theory, it is rarely calculated in this way because there are better ways of analysing it. The process of working out the shape of the IR from the input and output signals is called deconvolution and is generally more difficult. Therefore we are faced with a problem of deconvolution.
A very simple model of cloud feedback
Using this idea, one can construct a very simple model in which radiative flux heats the oceans, water vapour is generated and after some delay become clouds which modify the flux heating the oceans. This contains a number of simplifying assumptions, and is not intended to be an accurate model of the climate system but an illustration of how to go about using systems theory to model the climate system.
Assumptions:
1) There is no set point in this system. Control theory usually assumes some sort of demand and the system makes the output conform to that demand. This isn’t the case in climate.
2) The oceans are a perfect heat sink and, for simplicity they are well mixed within a very short time frame.
3) The cloud formation is described as a linear process; proportional to temperature and the clouds themselves have a linear effect on flux.
4) The heat required to produce a cloud in negligible compared to the heat in the oceans. This is a very important assumption from the systems theory point of view because the feedback, should not “load” the output.
One of the key words in these assumptions is linear. Most things in real life aren’t linear, but a widely used technique is to “linearise” a system. This assumes that for small change in an input, the output can be written as a series (Taylor’s) and that if the change is small, the output is proportional to the input because the other terms in the series are negligible. How far you can push this in any particular situation needs to be thought out carefully.
Given these assumptions, the model is consistent with the “thermodynamic” definition of feedback and the convolutions in the forwards and feedback paths can be analysed using a powerful technique called the Laplace transform, which is widely used in systems and control theory. Using this technique, the equation that describes the IR can be derived and is shown to be a second order differential equation (which from basic theory must be true).
|
Two quite interesting things come out of thinking about this very simple model. The first is that clouds could have both a positive and a negative feedback and the negative feedback has to be within certain limits to keep the system stable.
The second is that feedback systems may be prone to overshoot. A poorly damped system in response to a sudden step input will perform a series of exponentially damped sine wave oscillations above and below the level at which it eventually settles. Spectacular degrees of overshoot are unlikely, but small degrees are possible. It is important to recognise this effect because overshoot might be misinterpreted as “positive feedback”.
What does delay mean in a system?
Following this, how does one work out what is going on the system? A very important issue is delay. Although the static gain of a feedback system (i.e.: the climate sensitivity) is important, the dynamic behaviour of a system such as cloud feedback is determined by the patterns of delays imposed by the forwards and backwards components of the system. If there were a sine wave disturbance in a feedback system so that it was added to the output, negative feedback will cause a –ve sine wave to be summed at the input. However a negative sine wave is one that can either be thought of as phase shifted by 180 degrees or being delayed by half its period, so the negative feedback can be interpreted as a delay. If the delay were increased still further, the –ve sine wave feedback would become positive and the system would become unstable (or blow up). Suppose this system has perfect negative feedback at a particular frequency. What happens if we double the frequency of the disturbance? The feedback will be delayed by a whole cycle, rather than half a cycle, and the system is unstable. In practice, feedback is created by physical systems that store or dissipate energy, an in doing so, create a phase shift, or delay that depends on frequency. and so:
|
The impulse response of a system is related to phase shift because Laplace transform or (when calculable) the Fourier Transform of the IR is the system transfer function. This is a function of frequency, which describes how a signal at a particular frequency is attenuated and how it is shifted in phase. The Transfer function is an alternative definition of the system and is useful because it may be easier to calculate than the impulse response. Finally, convolution in the frequency domain is simply multiplication of the transform of the input signal by the transfer function and the output is the inverse transform of the result. This can be computed very efficiently.
Calculation of the phase spectrum of the impulse response of the model described above shows that it does not remotely approximate a pure time delay. While a time delay can be calculated for a particular signal, this will depend on the form of the input. Hence:
|
Both Spencer & Braswell, and Dessler use correlations to determine system parameters. A lagged correlation is simply a correlation between two signals, one of which has been shifted in time and is therefore a function of delay. This is a well-known signal processing operation, more formally expressed in the cross-correlation function, which is in theory incorporates all possible delays. An autocorrelation function is a signal correlated with itself, shifted by all possible delays. If one has an input signal passing through a system to produce an output signal:
|
It is not entirely clear to me what is being measured in the analysis by SB2011 and D2010, D2011 but when these methods are applied to the simple model discussed here, it is immediately clear that the parameters (gain and “feedback”) are a function of input signal bandwidth and so are not parameters at all.
How to measure the system parameters
One has to be very clear in what one is trying to achieve with an analysis of highly variable data and should start with a formal statistical hypothesis.
Hypothesis: The means of the distributions of the system feedback parameters are zero and the system cannot be distinguished from a system without feedback.
Having raised this hypothesis, one can explore the difficulties in testing it. A very common problem is that one can have an elegant mathematical model of a process but the data in too restricted to distinguish one’s model from a simpler one, in this case a model with feedback and the null hypothesis, one without feedback. This problem crops up in virtually every branch of science and generally results in an ill conditioned set of equations used to calculate the parameters of the model.
As a simple example, one has data of washout of a radioactive tracer in the body (appearing in breath or urine). This is of the form e-kt. One might say that actually there are two different processes causing the washout, one fast and one slow, so the decay curve could be 0.5.e-(k+a)t+ 0.5e-(k-a)t, ie.: two exponentials. These two models can be distinguished with “perfect” data but, with errors in sparsely sampled data, it is notoriously difficult to distinguish between them.
Another major problem is that if the model has an incorrectly modelled component, the real response of this component can “contaminate” the representation of another part of the system. In the model presented here, the ocean is represented as an integrator, simply because that is how SB2011 and D2011, D2010 choose to represent it. In practice, I bet it isn’t and it should really be represented as a multi-compartment system with flows and diffusion between different compartments. However, using real data to perform a deconvolution, the real behaviour of the ocean cannot be represented in the model and it will turn up spuriously in the “feedback”. This effect has been a huge problem in system identification in biology and medicine.
In climate because there is only one record of temperature, fluxes etc. and so one is dealing with essentially one observation, unlike classical systems analysis where one can “interrogate” a system using multiple input sequences and obtain a robust estimate of the system parameters. This may lead to formidable statistical problems in testing the hypothesis.
In practice, this depends on the length of the record, which determines the fundamental frequency of the discrete Fourier transform of the record and the length of the impulse response. In the time domain, this means that if the impulse response of the system is much longer than the record, one can’t characterise it. A further problem is that the Flux signals are band limited and may not fully cover the range of the impulse response.
Alternatively, if the impulse response is relative short, say a month, one will not be able to detect it with a record in monthly samples because it will be aliased[1]. ( For some reason, filtering a signal with a rectangular impulse response and then decimating the signal is commonly used in climate “to reduce noise”. This is highly undesirable because it may force aliasing on previously correctly sampled signal).
If we are lucky and the impulse response is relatively short compared with the length of the record, the record can be segmented into lengths that contain the impulse response and the parameters estimated from each record, allowing a better statistical estimate of the system response to be obtained.
An alternative approach, which I would use, is a parameter fitting technique. Although there are many schemes, the simplest is to use the model to create an output as function of the parameters and determine the integral squared error between the predicted and actual outputs. Minimisation of this integral with respect to the parameters of two competing models may determine whether the “best fit” model is one with or without feedback and allow the hypothesis to be tested.
The thesis I have advanced is that if one is to model a process, great care should be taken to ensure that the model actual conforms to the physical hypothesis, which should be stated explicitly and the assumptions underlying the model should be scrutinised. The data analysis methods should reflect the mathematics of the model and lead directly to estimation of its parameters. The use of multiple stages of analysis can become blind alleys and so one should ask what each processing step means in physical terms. I am somewhat sceptical that a systems approach of using deconvolution to imply feedback would yield an unambiguous answer, although it may be possible with better models of ocean temperature dynamics and cloud formation, coupled with a careful error analysis. I am, however, in no doubt that analysis of flux and temperature data, even at the modelling level presented here, is a formidable problem and requires a substantial amount of critical thought.
[1] Aliasing means that the sampling frequency is too low to capture the full frequency content of a signal. In this case, the signal is irretrievably corrupted and cannot be analysed. This is an error that has bitten many uncritical researchers in the hindquarters and given a sampled signal, the first question should be “is it aliased?”.
JC note: I have long been very interesting in getting perspectives from control theory experts on how climate science interprets the issue of feedback. I find this essay to be be very thought provoking and I look forward to your response. This is a technical thread and comments will be moderated for relevance.