Sunday, April 25, 2010

Plotting spatial trends in TempLS

This is a preview of Ver 2. It is a substantial extension, because it takes the principle of least squares fitting further, fitting spatial approximations to global patterns. Whereas before the spatial variation was attributed to local effects (station latitude, altitude etc), now the global model can be allowed spatial properties.

To do this, it is necessary to have a space of approximating functions. I'm using spherical harmonics - typically up to 200. So far I've analysed a model in which the global temperature is assumed to consists of a linear trend, spatially varying. But it will be possible to do yearly or monthly plots and, I think, trends and spatial averages for sub-regions, subject to enough data being available. I'd be interested in suggestions and comments, which is my excuse for posting at this fairly early stage.

Below the jump, I'll compare trends over the Earth's surface with plots from the AR4.

These plots used a total of 140 spherical harmonics, which gave reasonable resolution. The calculation is still quite fast, about 100 sec, so a lot more could be used. It is a whole sphere calculation - I've not plotted at the Arctic and Antarctic extremes, because the results then become extrapolations.

First, here are the plots of global trends from the AR4. It's Fig 3.9 from Sec 3.2.2.7.
and the caption says:
Figure 3.9. Linear trend of annual temperatures for 1901 to 2005 (left; °C per century) and 1979 to 2005 (right; °C per decade). Areas in grey have insufficient data to produce reliable trends. The minimum number of years needed to calculate a trend value is 66 years for 1901 to 2005 and 18 years for 1979 to 2005. An annual value is available if there are 10 valid monthly temperature anomaly values. The data set used was produced by NCDC from Smith and Reynolds (2005). Trends significant at the 5% level are indicated by white + marks.

The least squares method doesn't have those minimum requirements, although of course reliability diminishes when the data is sparse. Here are the TempLS plots for the global land/ocean fit:
I've tried to match the levels and colours, but it's obviously not perfect. The TempLS plots pick up the main features well, with an occasional extra flourish which may be real, since it can get information out of sparser data. One exercise will be to resolve this question.

1. Hi Nick,

I noticed in your post that you have a worldmap. Is that rworldmap available at CRAN R Project?

I've been trying for two days to make a plot with the legend that you have above and its driving me nuts.

I suspect it has to do with cutvector and the choice of colours. Could you put me out of my misery and show me how to code something like that?

SteveH :]

2. Steve,
I used the regular "world", which is the default in the map() function of the maps package. For the contours and key I used the filled.contours() function. Actually I hacked it a bit to get the legend down to the right size - I called my version fc(). Here's the code fragment I used to make a jpeg and a screen plot. G is the set of coefficients for the spherical harmonics that comes from the LS fit.

jet.colors <-
colorRampPalette(c("#00007F", "blue", "#007FFF", "cyan",
"#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000"))
ipcc.colors <-
colorRampPalette(c("#dddddd","white","#ffffcc","#ffddbb","#ffcc99",
"#ffaa77", "#ff9966", "#ff6666"));
ipcc.levels=c(-0.5,-0.2,0.2,0.5,0.8,1.1,1.4,1.7,2.4)*0.1;
#ipcc.levels=c(-0.8,-0.2,0.2,0.5,0.8,1.1,1.4,1.7,2.0)*0.5;

fillem=function(G){
lat=rep((lat0=-17.5:17.5 *5),each=72);
lon=rep((lon0=0:71 *5-177.5),36);
SK<<- array(P((p=spherharm(lat,lon, nSH)\$P),G), c(72,36));
SK[,c(1:7,35:36)]=0. # excluding the arctic
lon0=0:71 *5.3-157.5;
filcon=function(){fc(
lon0,lat0,SK,xlim=c(-160,lon0[72]),ylim=c(-180,180),
color=ipcc.colors, levels=ipcc.levels,
plot.axes=NULL,
key.axes=axis(4,ipcc.levels));
text(0,100,paste(c("Trend C/Dec",yr0+1,"to",yr1),collapse=" "));
};
jpg(paste("Trend",yr0+1));
filcon();
dev.off();
filcon();
}
I had to jiggle a bit to get the contour map and the world map to line up, hence the redefinition of lat and lon.

3. Expanding in spherical harmonics is a good idea. I assume you are allowing the spherical harmonic coefficients to vary over time, right?

That's how I've seen it done in other reconstructions (e.g., NRL-G2S).

4. Carrick,
This model just assumed the local seasonal etc behaviour, independent of year, plus a global linear trend. So the trend gradient varied over space but not time.
The next step is to model individual years (without the trend) with their own global SH coefficients. Probably not the whole time period - that would be a big calc, and unnecessary.

5. I guess the main issue I would have with that is that some regional patterns shift over time (e.g., ENSO).

6. Could you please do the saame for the period 1910 -1945 and for ocean end land sperate too?

7. 1910-1945 I could do, and I've realised that I should be comparing with GISS plots, which can be generated for any period. So there'll be a new post on that soon.

Separate ocean and land is trickier. The reason is that the spherical harmonic approximation method is inherently global. So if you do just ocean, say, it does provide a prediction for land as well. That doesn't matter too nuch, because it places no weight on deviation there.

8. fine thanks (sory for the typo's)

Ibrahim

9. I have the GISS file of station data after combining ( step 0)
so you can do an apples to apples

10. Thanks, Steven
I've picked that up. I'd need a matching inventory file (for lat/lon). There's one called v2.inv that came with the test run that Ruedy provided for CCC - looks like it might match.