Guest Post by Willis Eschenbach
I’ve been thinking about thermal lags in the climate system. Everyone is familiar with thermal lag in everyday life. When you put a cast-iron pan on the stove flame, it doesn’t heat up instantly. Instead, the warming process starts instantly, but it takes some amount of time to actually up to a stable heat. And when you turn the heat off, it doesn’t cool down immediately. That’s thermal lag.
Now, Dr. Shaviv claimed that I didn’t understand harmonic analysis and the lag between the sun and what it is heating. I’ll return to that question further down. For now, thermal lag is generally modeled as an “exponential” lag. In an exponential lag, say when you take a hot pan off of the flame, it doesn’t become cold immediately. At first it cools fast. Then as time goes on and it nears room temperature, it cools slower and slower. That “fast at first, slower later” form is characteristic of an exponential process.
In the calculation of an exponential lag, the delay in the response is governed by a time constant called “tau”. The larger the time constant tau, the greater the lag time, and the slower the pan will either heat up or cool down.
With that as prologue, here’s a graph of the heat being turned on and then off under the pan on the stove (black/yellow line), and the temperature response of the pan itself (lines colored blue to red).
Figure 1. Example graph showing how a pan heats up and cools down. Time is in minutes. The fire is shown by the black/yellow line. The fire goes on under the pan at about time t = 8 minutes, and comes off again five minutes later at about t = 13 minutes or so. The time constant “tau” is measured in whatever the underlying time unit is, in this case minutes. Amplitude in arbitrary units.
Figure 1 shows graphically the change in pan temperature, using different colored lines for different values of the time constant tau. The blue lines show small values of tau. A small value for tau would be the situation with say a small thin aluminum frypan, which would rapidly come up to heat and stay at temperature until the fire is turned off. At that point, it cools down quite rapidly. We say that it has small “thermal mass”, so it heats up quickly.
At the other end of the spectrum, the red line shows a large time constant tau of 9 minutes. It would be characteristic of say a big, heavy cast iron cook pot. It would have lots of thermal mass. In the example above, note that the iron pot heats up so slowly that it doesn’t even get up to heat before the fire is turned off.
This illustrates at least three points about this kind of exponential thermal behavior—exponential action, changing shape, and changing size.
First, the action happens exponentially. This means quickly at the start, and then slowing down as time goes on. For example, a hot pan taken off the fire cools rapidly to begin with, but as it approaches room temperature, the cooling slows.
Second, and this is an important point. In general, the shape of the exponential lagged response (the pan temperature) is NOT the same as the shape of the impulse (the on/off of the fire). In the red lines in Figure 1 above, you can see how the square-wave on/off thermal impulse of the flame (black/yellow line) is transformed into a kind of shark-fin shaped thermal response in a heavy cast-iron cookpot.
Third, if the time of the impulse is short enough, the amplitude of the response is smaller than that of the impulse. In other words, in Figure 1 the blue line gets as high as the black/yellow line, but the red line doesn’t get that high before the time comes when it starts dropping again.
There is a special case, however, where the exponentially lagged response has the same shape as the impulse. This is when the impulse has a sinusoidal (sine wave) shape. An example of this would be the solar strength as the earth moves nearer to and further from the sun over the course of the year. Over the course of a year, solar strength varies sinusoidally about 23 W/m2 on a 24/7 global average. It peaks in early January, and bottoms out in early July.
Anyhow, I got to thinking once again about the effect of exponential lag on a sinusoidal input. So I made up a graph like Figure 1, except using a sinusoidal input instead of a square wave input.
Figure 2. Lagged versions of a sinusoidal input. The input is the black/yellow line in the back of the drawing, which represents a flame under the pot of varying strength. The lagged versions showing temperature variations in the pan are shown in color. Time units are quarter cycles. If the full cycle is one year, tau is in months.
Now, I was very happy when I first saw that result. I was elated because I’d never been able to calculate the thermal lag directly from the observations. Oh, I could fit a lagged curve to the results, but to me that’s unsatisfying. And in this graph, I noticed a very curious thing that I’d not known before.
This is that the peaks of the lagged versions of the input wave fall exactly on the original sinusoidal impulse line (black/yellow). This was important to me, because I’d never been able to figure out how much smaller the lagged response wave was than the size of the original impulse. For example, the peak in the earth’s temperature lags behind the peak in the sun. This makes the thermal response smaller … but how much smaller? The graph gave me hope that I could calculate it, and after much investigation of the various questions, much walking in the forest and lots of headscratching, I’ve figured out a heuristic formula that relates all three variables—the time constant tau, the observable lag in the timing of the peak temperature versus the peak of the input, and the scale factor reflecting the loss of amplitude due to the lag. I’m sure some better mathematician than I can provide the actual functions, but the ones I have are very accurate. I’ve appended the functions at the end.
And there is a lot more to learn from Figure 2.
One thing is that the longer the lag, the smaller the resulting thermal response. Expected, but good to see just how it works.
Another is that the resulting lagged signal is also a sine wave. In this special case, the exponential lag doesn’t change the underlying sinusoidal shape of the input impulse (black/yellow line). An increasing time constant tau just pushes the sine wave further and further back in time after the impulse, and scales that sinusoidal impulse ever smaller with increasing lag.
This has an important corollary. If we are looking for the result of sinusoidally varying forcing, the thermal response will have the same shape as the input forcing, but it will occur later in time.
An important observation is that with increasing tau (blue to red, tau of 1 to 9), the colored lines get closer and closer together. As a result, no matter how big the time constant tau gets, the observable lag will never be greater than a quarter of a cycle. For a yearly varying input like the solar energy, that means that the peak temperatures in the ocean and the land can never lag more than three months behind the peak solar … and if they did they’d have zero amplitude. Figure 3 shows the lags for values of the time constant from one to nine.
Figure 3. Closeup view of the upper left section of Figure 2. Vertical black lines show the time lags in the peak of the observable response with values of tau varying from one to nine. Time units are quarter cycles. If the full cycle is one year, tau is in months.
Notice how the vertical lines bunch up at the right. No matter how large tau gets, the actual lag in the response peak is always less than a quarter cycle.
So now I have have developed the tools do the lag math. I can take a look at the lag between solar input and thermal response, and from that I should be able tell you the time constant tau, as well as the scale factor lambda, which relates the solar input to the change in temperature. Let me put this to a bit of practical use.
I looked for starters at the Northern Hemisphere ocean temperatures, and the corresponding NH sun. I’ll use the data from the CERES satellites. The data is monthly. To increase the accuracy of the calculation of the observable lag, I first fit a spline to both monthly datasets. I spline them at 360 points per year, then I average the points. I do this for both the solar data, and for the ocean temperatures. Then I graph them against each other as follows:
Figure 4. Scatterplot of the monthly average northern hemisphere solar anomaly versus the monthly average of the northern hemisphere ocean anomaly . Data has been splined before averaging with a resolution of 360 points per year.
So here’s how the procedure works. Using the splined version shown in Figure 4 allows us to measure the lag between the solar forcing, which has extremes in December and June, and the lagging ocean temperatures with peaks a couple months later, in February and August. This is about two months. Splining it to 360 points puts the measurements in degrees. There is a 61° lag in the summer, and 58 degrees in the winter. So the lag is almost exactly sixty degrees.
Then, using the relationships established above, I find the following:
Observed lag angle: 61° summer, 58° winter
Time constant tau: 3.3 months
Tau-dependent lag scale factor (amplitude reduction from no-lag condition); 0.48
Thermal sensitivity lambda factor (accounting for lag scale factor); 0.067°C per W/m2
I have not mentioned the “lambda” factor in the lagging equations. This is the transfer function or value that converts from the range of solar energy input swings to the range of the resulting temperature swings. Lambda in this case is .067°C per watt per square metre (°C per W/m2). This means that when the solar energy changes by 100 W/m2, the ocean changes temperature by about 7°. Note that this thermal sensitivity (7 degrees per 100 W/m2) is about twice that shown in Figure 4 above. This is because of the lag scale factor—the lagging alone cuts the response about in half, which means that the actual sensitivity is about twice what the observations show.
Now, having calculated those values, I used them to do my test of whether my conversion functions were correct. These were my functions to convert between observed lag, tau, and lag scale factor. Holding my breath and hoping I’d squashed all of the bugs, here’s what I got:
As you can see, it appears that my method is quite accurate. The values for the time constant tau and the lag scaling factor were calculated directly from the one observation, that of the observed lag between the peaks. It was done without using a fitting process of any kind. This indicates that my functions are, well, functioning well …
Next, let me compare the thermal response of the northern hemisphere land with that of the northern hemisphere ocean.
Both the land and the ocean get the same sunlight. As you would expect, the land temperature varies much more than the ocean temperature. Here are the calculated thermal lag parameters for the land, including ocean comparisons.
Observed lag angle: 25° summer, 22° winter (ocean 61° summer, 58° winter)
Time constant tau: 1.0 months (ocean 3.3 months)
Lag scale factor (lagged amplitude divided by no-lag amplitude); 0.91 (ocean 0.48)
Thermal sensitivity factor lambda (accounting for lag scale factor); 0.137 (ocean 0.067°C per W/m2)
Why is the thermal sensitivity factor lambda higher for the land than for the ocean? The answer is, they vary in the thermal mass which is affected by the yearly swings. Curiously, this can be calculated directly from lambda, the thermal transfer function. In a twelve-month dataset, where tau is measured in months, the mass involved in the swings is the mass of a depth of ocean in metres which is 8 times tau. In other words, on land tau = 1 so an amount of mass equal to about 7 metres of ocean is involved in the annual swings. In the ocean, tau = 3.3, meaning the water down to 3.3 * 8 = 26 metres is involved in the annual swings.
Now, I speculated a while back that the lag at the daily level is related to the lag at the annual level. I just realized that I can now answer this question by simply examining a graph of the hourly surface temperature where I live. If my speculation is correct, the time constant for the individual station data should be on the order of that for the entire land mass. That was a tau of 1.0, which is 0.52 radians. I wouldn’t expect too great a fit, because one is a point measurement and the other is a hemispheric average. Anyhow, here’s the Santa Rosa data.
Now, using the tools I’ve developed I can take just the information visible in this graph, no digitizing, and calculate the time constant tau for that situation. Tau in Santa Rosa for the daily data turns out to be 0.48 radians, compared with 0.52 radians for the monthly values for the entire NH land …
So it does appear that the time constant (in radians) is period-invariant, meaning that it is related to the lag angle and it doesn’t vary with the length of the cycle.
• It is possible to use solely the observed lag between an impulse and the thermally lagged reaction to directly calculate the time constant tau as well as the scale factor resulting from the lagging alone. I’m sure there are mathematicians out there who both knew that and can derive the formulas mathematically.
• Knowing tau and the lag scale factor allows us to calculate the thermal transfer factor lambda.
• In exponential decay the time constant tau is also known as the “e-folding time”. This is similar to the concept of “half-life”. Half-life is the time it takes for something to decay to 0.5 (half) of the original value. Similarly, e-folding time is the time it takes for something to decay to 0.37 (1/e) of the original value.
Let’s see what that means in practice. Over the land, the time constant tau is about one month. That means that after a month, the original signal is down to 0.37 of the original value. So far so good.
So where is the signal after two months? Well, another month has passed, so it is now 0.37 of its value at the start of the month, which is 0.37 times 0.37 = 0.37^2. And so on down the line. This means that it is dropping off fast. In general it falls off as e^(-t/tau). In practice, this means that with a time factor tau which is 2 * pi / 3 or less (a tau of 4 in a twelve-month cycle), any signal will decay to within a percent or two of zero within one cycle.
Since I’m looking at observational data for the temperatures, the data contains all possible responses to the rhythmic solar forcings. These include but aren’t limited to resonances, reactions, changes in clouds, feedbacks, slow and fast exponential increase and decay from past inputs, regime shifts, variations in water vapor and thus the overall greenhouse effect, the speed of the throughput of the energy, and all of the rest.
Now, what are we measuring with the thermal transfer coefficient lambda? We are measuring what the thermal response of the system would be if there were no thermal lag. But if there were no thermal lag … then it must be an equilibrium thermal transfer coefficient, or quite close to one. This is interesting because this implies that the lambda values for the no-thermal-lag condition must be close to what is called the “equilibrium climate sensitivity”. This is usually expressed as the temperature variation arising from a forcing change of 3.7 W/m2, as this is the calculated change in forcing from a doubling of CO2. For the ocean this gives us a “climate sensitivity” of .07 °C per W/m2 * 3.7 W/m2 = 0.26 degrees per doubling. The land figures are 0.14 °C per W/m2 * 3.7 W/m2 = 0.5 degrees from the forcing change from doubling.
Now, I see this as a measure of the strength of the climate thermal regulation system. The average swing in incoming solar energy in the Northern Hemisphere is two hundred watts per square metre peak to peak. The average swing (weighted ocean/land average) in no-thermal-lag NH temperature is 7°C/.48 lag factor * .7 + 26°C/ .9 lag factor * .3 = 19°C peak-peak. On the other hand, If we figure the blackbody temperature change from that 200 W/m2, if there were no thermal lag that would swing the northern hemisphere’s temperature by 37°C peak to peak … and instead it only swings about half of that, about 19°C. This is from a combination of fixed and variable ground albedo, variable cloud albedo, emergent thermoregulatory phenomena, biosphere variations, and all of the other ways that the climate responds to the changes in forcing. And it’s a good thing we do.
Finally, is there a long-term residual that keeps going for a long time? Possible, to be sure. For example, there could be 11-year cycles from sunspots, or longer cycles from the swings in solar activity. Dr. Shaviv is looking for such longer cycles in the satellite sea level data, and I wish him well. It’s a crazy world out there, and I’ve been surprised before at what people discover. However, given the small size of the thermal time constants, in the range of one to four, it seems such signals would die out for practical purposes in fairly short order
My best regards to everyone, 2:17 am of a foggy night. On the radio, Roger Miller playing through the headphones, “I’m a man of means, by no means … king of the road” … it’s that kind of night for me, and I hope for you as well.
[UPDATE] the always reliable Joe Born writes in the comments:
<blockquote>When you’re talking about a one-pole (“one-box”) system, as Mr. Eschenbach implicitly is (the climate systems isn’t one, but for his purposes we’re assuming it approximates one), the attenuation theoretically is 1/sqrt(1+(tau * omega)^2), and the lag is theoretically arctan(omega * tau)</blockquote>
I’ve used these to re-write my formulas. My old heuristic formula gave the tau of the ocean as 3.27 months … the actual formulas from Joe give a value of 3.30 for the tau of the ocean. So my conclusions are still solid.
I replied to Joe:
Joe, here’re my new formulas based on your equations. They all take native units (e.g. hours, months) for the lag (peak forcing to peak response), tau, and the period. So a function call of lag2tau(2,12) means 2 months lag in a 12-month year, with the result being tau in months.
Attenuation is a decimal from 0 to 1. The empty function calls call the function with the default values, it’s how I test my functions.tau2lag=function(tau=3.3,period=12) atan(tau*2*pi/period)/(2*pi/period) tau2lag()  1.998003 lag2tau=function(lag=2,period=12) period*tan((2*lag*pi)/period)/(2*pi) lag2tau()  3.307973 tau2atten=function(tau=3.3,period=12) 1/sqrt(1+(tau*2*pi/period)^2) tau2atten()  0.5009052 lag2atten=function(lag=2,period=12) tau2atten(lag2tau(lag,period),period) lag2atten()  0.5 tauconvert=function(tau=3.3,period=12,newperiod=24) tau/period*newperiod tauconvert()  6.6
Note that omega is equal to 2*pi/period, again with the period in native units. Also, since tau is expressed in units of omega it is not conserved when you go between say a 12-month and a 24-hour situation. Use the “tauconvert” function to convert tau between periods. I’ve deleted the old R functions at the end.
My Perennial Request: If you disagree with someone, please have the courtesy to quote the exact words you disagree with. This allows everyone to understand exactly what you object to.
Math Notes: The main equation I use calculates the change in temperature per unit time as
where the subscript t shows time, T is temperature, F is forcing, tau is the time constant, and lambda is the transfer function or coefficient.
Adopting the terminology that for any variable X, the variable ∆X is defined as:
the main equation simplifies to
The first term on the right hand side, involving the forcing change ∆F, is the direct action of the forcing. “Lambda” is the transfer function, which scales the input to the output. The second term on the right, involving the temperature change ∆T, is the lagged action. Incoming forcing is partitioned between the two branches (immediate action and delayed exponential decay action) by the alpha and (1-alpha) terms. As tau gets larger, the immediate direct action decreases and the length of the lagged action increases.
If tau is zero, alpha is zero, and there is no lag. The response is simply the input times lambda, the transfer factor that relates the range of the thermal response to the range of the input forcing.
Tau and lambda are the only parameters in Equation 1. My goal in this investigation was to be able to calculate tau and lambda directly from nothing more than the delay between impulse and response. In that, I’ve been successful.
Dr. Shaviv’s Explanation
I said above that I’d return to Dr. Shaviv’s objections to my previous postings, which I append below:
Eschenbach has no clue on what the meaning of harmonic analysis is. He still thinks I assumed the sea level should be in sync with the forcing…
Suppose the solar forcing is Delta F = F0 exp (i 2 pi t/P)
The heat in the oceans (per unit area) will be (assuming feedback is small)
Q = int Delta F dt = – i P F0 / (2 pi) exp (i 2 pi t / P)
The sea level, if it is only affected by thermal expansion, will then be
Delta h = alpha Q = – i alpha P F0/(2pi) exp (i 2 pi t / P)
Delta h = alpha Q = – i alpha P/(2pi) Delta F
It seams like an algebraic relation, but in fact, the -i gives a 90° phase shift which describes the integral(!) If one take the real part of the above equation (to get the physical quantities), one finds
Delta F = F0 cos (2 pi t/P)
Delta h = alpha P F0/(2 pi) sin (2 pi t / P ) = alpha P F0/(2 pi) cos (2 pi t / P – pi/2 )
That is, the integral appears in a harmonic analysis as a phase lag of 90°. In real life, this angle is not exactly 90° because of additional processes taking place, which was the whole point of the analysis!
So again, Eschenbach’s analysis, compared the sea level to the solar forcing and found there is a poor fit. It is because the sea level is supposed to lag by 90° after the solar forcing! Eschenbach still does not acknowledge this.
Also, Eschenbach’s claim that I compared to the sea level and not the rate ignores the fact that a harmonic analysis is expected to give a 90° lag if the sea level is the integral of the forcing, and this is exactly the case.
OK. I’ve corrected a couple of Dr. Shaviv’s math typos, in red. I hate math typos, easiest thing to do, and I certainly hope I don’t have any in my math above.
There are two problems in the above explanation. First, Dr. Shaviv is discussing two kinds of lags, without noting the distinction. Second, except for in the math, he uses “forcing” to mean “CHANGE in forcing”.
I’ll start with the question of the two lags, and I’ll use his terminology from his math. The first lag is the 90° lag in a sinusoidal system between the change in forcing Delta F (∆F) and the thermal response Q. The second lag is the thermal lag discussed above. The first sinusoidal lag is a result of the curious fact that the integral of a cosine is a sine, and a cosine is just a sine shifted by 90° … a couple facts that I’ve known since I taught myself calculus in high school.
Note that this lag he is discussing is between ∆F, the change in forcing, and Q, the result. However, in a sinusoidal system the same lag exists between ∆F, the change in forcing, and F, the forcing itself. In other words, although there is a 90° lag between ∆F and Q, there is the same delay between ∆F and F. As a result, there is no lag at all between F and Q.
Now, recall that Dr. Shaviv said above that:
So again, Eschenbach’s analysis compared the sea level to the solar forcing and found there is a poor fit. It is because the sea level is supposed to lag by 90° after the solar forcing!
Here you see the second problem. Dr. Shaviv is incorrectly using the term “forcing” F to mean the CHANGE in forcing ∆F, the quantity which in his equations above he calls “Delta F”. Yes, sea level is supposed to lag Delta F by 90°. But if there were no thermal lag, there would indeed be no lag between sea level and forcing F.
This is a continuing problem with Dr. Shaviv’s claims, which I’ve pointed to before. For example, he also says above:
… a harmonic analysis is expected to give a 90° lag if the sea level is the integral of the forcing, and this is exactly the case.
Not true in the slightest. Look at his own math. In his analysis above, he doesn’t use the integral of the forcing (which would be “int F dt”). His exact statement was:
Q = int Delta F dt
He says that the change in sea level is alpha Q, and he says this sea level change 90° out of sync with the solar forcing F.
But take another look at his equation immediately above … what is the integral of Delta F? … I didn’t see it the first time through myself, but as you might expect the integral of the derivative of the forcing F, which is his int Delta F dt above, is the forcing F itself. That means that Q = F, and it means that in the absence of thermal lags, the sea level changes and the forcing indeed are in sync.
So when Dr. Shaviv’s says that “sea level is the integral of the forcing”, he is simply wrong. His own math shows that sea level is NOT the integral of the forcing F. It is the integral of the change in forcing Delta F ∆F.
The other kind of lag is the lag I discuss above, the lag between the peak forcing F and the peak thermal response. If there is no thermal mass then the thermal response is immediate and the lag is zero. As shown above, the thermal lag can never be more than a quarter cycle, and the thermal lag is a function of a time constant tau.
With those two kinds of lags in mind, consider Dr. Shaviv’s statement:
That is, the integral appears in a harmonic analysis as a phase lag of 90°. In real life, this angle is not exactly 90° because of additional processes taking place, which was the whole point of the analysis!
The phase lag between cos(x) and int (cos(x) is 90°. Always. Exactly. In real life, it is still exactly 90°. It is not changed by “additional processes taking place”. It is a mathematical fact. However, it is not much considered, because in real life the forcings and the responses are NOT 90° out of phase as Dr. Shaviv claims. They are generally in phase, with the response slightly lagging the forcing. If the sun goes stronger and weaker, and the sea level is a function of the sun, we expect the change in sea level to be somewhat thermally lagged behind the forcing … not 90° out of phase plus or minus something.
Let me be clear that I think that Dr. Shaviv knows the difference between the forcing F and the change in forcing ∆F. Otherwise, he couldn’t do the math. I don’t know why he keeps confusing them over and over in his responses.
Finally, I was amused by a commenter on Dr. Shaviv’s quoted passage who said:
I’d rather [Willis] took some math classes. His previous teacher was not up to the task. (He identifies as self-taught). He’s just a few semesters away from a “who knew” moment.
It’s true, I’m self-taught in math. When I was a junior in high school, I finished all the math courses that the school offered. So my math teacher, Mr. Hedgji, bless him, gave me a calculus textbook and said “Go for it”. He helped me when I had questions, and otherwise left me alone to go for it.
So I did, and I never stopped.
I notice that the anonymous commenter couldn’t find any actual math errors of mine to comment on, so he resorted to personal judgements. Meanwhile, I’m the one who puts my math work out here on a regular basis so that folks like the charming commenter can snipe at it, and so that serious mathematicians can point out my flaws. And I’m the one discussing the problems with Dr. Shaviv’s math. Not only that, but I give both the commenter and Dr. Shaviv full opportunity to find any problems with my work and point them out in public. Strange as it may seem, I welcome people finding errors in my work. Folks like Mosher and Nick Stokes and Dan Hughes and too many other commenters to name have saved me literally years of wasted effort running blind alleys by pointing out logical, mathematical, theoretical, or other errors in what I was doing, and my hat’s off to them for their due diligence. It’s no fun at all to have my errors pointed out in public … but it’s far better than the alternative.