Just a short note today on code for McShane and Wyner and its comments. The discussion of McShane and Wyner has been greatly enhanced first by their archiving of code and secondly by simplifications of the code both by the online community and discussants.
Most of the McShane and Wyner analysis is in R, making this part of their analysis quickly accessible to a wide community. (Their Bayesian analysis is not in R – it seems to be in WinBUGS.) There’s a lot of self-teaching in R programming and I nearly always learn something from reading other people’s code. Sometimes you see things that you think that the other guy could learn. For example, I find it handy to organize R data sets as time series (“ts”) objects – not used by Mc-W. I like to have objects in use available on the internet for download using download.file(…, mode=”wb”), rather than packing things up in zip-files. It’s a small stylistic preference.
I was able to get the R portion of their code to operate quickly and easily, but I’m not familiar with WinBUGS and wasn’t sure how much overhead would be involved in figuring it out.
Meanwhile, I spent some time looking at the code for Gavin Schmidt’s (Schmidt, Mann Rutherford) comment on McShane Wyner. The code is night and day different from the code of Mann et al 2008, which is ponderous and opaque to the point of virtual unintelligibility. While the entire point of the Schmidt comment is nothing more than preserving bristlecone impact without mentioning bristlecones, the code of their comment is clean and easy to follow.
One useful improvement to the McShane-Wyner code is the collection of their code into a well-organized master script, much improving the ease of replicating results. Their master script makes extensive use of the source command. Occasionally they have second-tier source commands (called from a source script) rather than keeping all the source commands in the master script, but this is a stylistic nit.
Their other major improvement is applying the port of the Bayes portion of McShane-Wyner to R originally developed at a blog blog here. This port requires the R-package rjags AND the separate installation of JAGS 2.1.0 . It took me a while to figure out that JAGS 2.1 needed to be installed separately, but once that was done, I was able to proceed uneventfully.
For the most part, the reconstructions illustrated in Mc-W Figure 14 are based on a variation of MBH98-99 methodology (and not the lasso methodology used in the red noise analyses) i.e. they did a principal components analysis of proxies and gridded temperatures (with the M)8 network) and reported various permutations, yielding a lot of different squiggles. The big issue in the various permutations is not discussed by Mc-W or (especially) by Schmidt, as it pertains to the extremely stale issue of the presence/absence of strip bark bristlecones.
I mentioned previously that the 93-series network used by McShane-Wyner included 24 strip bark chronologies and 2 Korttajarvi sediment series. Because there are so many strip bark chronologies and no attempt is made to downweight all these nearby chronologies, they dominate the PC1 in the AD1000 network even without Mannian centering – a point originally made in McIntyre and McKitrick 2005 (EE), but mostly ignored. It can therefore hardly be a surprise that the “red” hockey stick shaped reconstruction in McW Figure 14 (see below) looks like MBH99, because it is simply one more avatar of MBH99 bristlecones.
In the following graphic, I’ve indicated the weights of individual gridcells by the weights of the filled circle – a technique used here on a number of occasions. After all the huffing and puffing of megaflops of multivariate analysis, it all comes down to a vector of linear weights of the original proxies (in calculating the weights, I’ve used linear algebra shown at CA in the past.) Other permutations are less weighted towards bristlecones.
Schmidt et al contested the Mc-W 93-series network as not being “properly” selected. They reduce the network to 55 proxies, excluding the two Tiljander series retained in the MW 93 and tree ring series with fewer than 8 cores. Bristlecones remain well represented, needless to say. I’ll discuss their version on another occasion, together with a no-bristlecone version.
The availability of code stimulates discussion since you can see what they actually did and it hugely mitigates the arid controversies all too typical of climate science about whether one had correctly guessed at some bizarre methodological quirk.