As I discussed previously here, I had been able to quite accurately replicate Mann’s temperature principal components (the first 16 of which have been archived for a long time). Notwithstanding that these calculations can actually be replicated quite accurately, Wahl and Ammann calculated new reference temperature PCs using annual data on the basis that they get similar results. Perhaps. But if the purpose is replication, it makes a lot more sense to use the MBH98 temperature PCs (or at least replications of them) . Also pardon me for being a bit skeptical about Hockey Team assurances that something “doesn’t matter”. I’m getting to find my way around WA code a bit better now and I substituted the Mann PCs for the WA PCs and carried on. The normalized proxies used by WA were very slightly different than mine, owing to slight differences in the order of rounding and normalization: in some cases, they rounded before normalizing, while I didn’t (I think that my interpretation of MBH methods is the more probable).
Then I re-did the calibration-estimation calculation up to Reconstructed Principal Components. In Early Returns, based on Wahl-Amman PCs, I reported that there was a correlation >0.999 between the two calculations, but scaling differences. Well, using MBH98 PCs and non-rounded proxies, guess what: the replication of the RPC1 was exact. I mean to 10 decimal places for all 581 values!
You’d think that this would have been worth a mention in their press release and article. I’m still wading through the downstream scaling to NH temperature. Because the scales were all different using the WA PCs, it complicated this process considerably and I didn’t finish that yesterday. However, now that I’ve got everything on the same scales, de-coding the NH scaling should come a little faster today. For reference, the incisions were surgical on the following: svdproduct$TmatU, svdproduct$TmatV, svdproduct$TmatD (temperature PCs) and datamatrices$Ymat.anomaly (proxies).
I’m working on the 15th century step only for now. I’ve unpacked the functions into a running script.
Afternoon Update: I think that I’ve pretty much reconciled the scaling to NH temperature.
This has always been a bit mysterious – Ammann’s code shows a scaling step added in April 2005, so he seems to have had some trouble here as well. Von Storch et al. seemed unaware that there was a scaling procedure in MBH98 – I mentioned this in an earlier post [here]. In our previous calculations, because I couldn’t get the scaling to reconcile, I ended up simply re-scaling the emulation to NH temperature. I should be able to eliminate this plug in my emulation, which will be nice. I anticipate being pretty much compatible with WA by tomorrow. I’ve done some preliminary runs and obtained R2 and RE results pretty much identical to previously reported results. Once I’ve got complete replication, I’ll re-execute all reported results and report on them. That may be a few days as I’ve got a few other things to do simultaneously. It’s interesting to see some programming approaches in WA. I like some ways that they handle output logfiles.
On the other hand, their handling of matrices is often very awkward. For example, there is a profusion of differently named matrices such as “Tmat.anomaly.cal.fitted.NH.mean.corrected” and “Tmat.anomaly.precal.fitted.NH.mean.corrected”. The phrase *.cal.* and *.precal.* distinguish matrices in the calibration period 1901-1980 and the “pre-calibration” period 1400-1901. These could have been combined into one matrix using an index for the calibration period (or other period) e.g. Tmat.anomaly..fitted.NH.mean.corrected [cal,]. This would eliminate needless duplication. I also don’t like their handling of situations where each column of a matrix is multiplied by an element of vector. For example, if A is a matrix and v is a vector, the simplest way of doing this is A %*% diag(v), (%*% is matrix multiplication). WA typically will do t( t(A)*v). It works but you lose sight of what’s going on, especially when they create new and unnecessary matrix entities with *.transpose suffix. These are not big criticisms, but they are matters that could have been done a little simpler. (I’m sure that there are lots of criticisms of my coding style, but these are just things that I noticed and had to think about to figure out what they were doing.)