Adam B. Smith

Importing NLDAS and GLDAS data into R

August 29th, 2016

NLDAS/GLDAS

OK, I just spent the entire day obtaining and learning how to import the NLDAS and GLDAS data into R. This could have been made a lot simpler with better meta-data descriptions and “readme” files placed in locations they need placed. In any case, I’m posting this to save anyone else wanting to use this data some precious time. In case you didn’t know (I didn’t until yesterday), the NASA Land Data Assimilation Systems (NLDAS) and Global Land Data Assimilation Systems (GLDAS) are measured/interpolated and/or modeled climate and land surface variables for essentially the conterminous US (NLDAS) at 0.125 deg resolution or the world (GLDAS) at 1 deg resolution. There are a lot of variables of interest, including the basic set of min/max air temperature and precipitation, plus snowfall, soil temperature, LAI, albedo, incoming/net shortwave and longwave radiation, etc. The data is available in sets representing calculations at 3-hr intervals or monthly intervals or averages for each month across the given time period. Most of the temporal extents of these models cover 1979 to the present.

Both NLDAS and GLDAS have version 1 and 2, the latter being newer and more sophisticated. Both have also been run with 3 land surface models: Mosaic, Noah, VIC–but wait, there’s a fourth, SAC, which is not described in the ReadMe file for NLDAS2. There are also “FORA” and “FORB” data sets put alongside the three models with little explanation as to what they are. They contains the forcing variables used by the land surface models. Note that the forcing variables remain unchanged across the three land surface models, so if you want a variable in the forcing set you can just get it from FORA or FORB.

Obtaining the data

If you want the entire hourly dataset it will take a long time to download as each model set contains tens of thousands of files. The monthly is only a few thousand; the monthly averages only 24 (one raster file plus one XML file per month)–these have the word “climatology” in their data set names.

There are many ways to get the data, including using wget which snatches files from a list of links but whose help is written for Unix/Linux, which is a girl I used to know. You can also get subsets of the hourly data using the Simple Subset Wizard (search for “NLDAS” or “GLDAS”). The SSW can also export the files in NetCDF format, which obviates the stuff below but I found the SSW did not always give me the full set of results. Hourly/monthly data are available from NASA’s Mirador (same search) or GES DISC. I used the latter then DownloadThemAll, a Firefox plugin. This still took a lot of clicking, but not near as much as if I had done it manually. (You’ve also got this badly documented FTP site.)

Extracting the data

So… the problem is that G/NLDAS files are stored in GRIB format, akin to a raster brick, but with no meta data on layer identity that is automatically imported into R or ArcMap when the raster is loaded. The XML file that comes with each GRIB file has a list of variables, but they are not in the order they appear when imported into R. Likewise, when imported into R the metadata that should come with a GRIB file is not associated with the file contents, so you are left with a long series of rasters with many meaningless numbers.

Solution:

In R:

library(rgdal)
library(raster)
grib <- readGDAL(‘<filename of GRIB file>’) # read GRIB file
grib <- brick(grib) # convert to raster brick
grib # notice brick has N layers

2. Download wgrib.

3. Open a command (DOS) window and navigate to the folder with wgrib. In Windows you can get a DOS window by pressing the Windows key then typing “cmd”.

4. Issue “wgrib <filename with no spaces>”. The output will show a table with variable names and attrbutes for each layer. You will need to copy the GRIB file into the same folder as wgrib or put it into a folder with no spaces in its name or any of its parent folders. Probably a way around this…

5. Consult the metadata file “README.NLDAS2.pdf” from the G/NLDAS website and see Table 4a therein. Find the “Short Name” of the variable you want.

# 6. Now look for that variable name in the DOS command window. The output from wgrib will show you the layer number of that variable. Remember this number… call it “x”.

7. Back in R:

myLayer <- grib[[x]] # the layer you want

NB This seems to work for every variable except TSOIL (soil temperature) which for the file I experimented with has 3 such layers. I am guessing these pertain to the three soil layers for the particular land surface model. There was also a band named “var255” at the end, which had what seem like meaningful values of some variable.

Note that the layer you want may not be in the same place across land surface models–i.e., LAI may be layer x in one and y in another.

As my professor said to me once, if things don’t go well for you at least make it better for the next person.

Posted in Nichey, NLDAS and GLDAS | No Comments »

Weighing the importance of scale

August 1st, 2016

I just finished an exciting read: Schweiger & Beierkuhnlein’s study on how well temperature predicts distribution of 19 vascular plants across 3 spatial scales (ranging from ~<1 m to 1000s of km). Overall they find regardless of the scale the same optimum temperature is observed (weak scale dependence). Nonetheless, they also find that the maximum probability of occurrence increases with grain size (strong scale dependence), which they interpret to mean that temperature is a more important driver of distribution at coarse scales.

Cool stuff! I’ve been specifically wondering about this for a while, and this seems to be the first test thereof.

But what does it really mean? First, I should say that their analysis was based on extracting metrics from the modeled response curves, not the response curves per se–and to my eye the curves for any particular species seem very different across scales even after correcting for differences in height (their Figs. 1 and S1). I would have liked to see a statistical comparison of the shapes of curves.

But let’s let that lie and think about what they found. In a nutshell, their results are predicted by the Eltonian noise hypothesis which posits that abiotic drivers like temperature will be more important at coarse scales while biotic drivers will create “noise” in distribution at fine scales–noise that will be generally imperceptible at coarse scales. They infer this from the fact that maximum probability of presence increases with coarseness of grain (i.e., when predicting presence at fine grains the maximum probability will be low). Ergo, temperature is a stronger predictor of presence at coarse scales.

While I can’t refute this observation on face value, I do wonder if the maximum probability of occurrence that they estimated at coarse grains is more than expected by chance based on combining probabilities of presence at fine grains. Consider for example, a simple situation where the “coarse” spatial domain (of area A) is composed of 2 fine-grain domains (each of area A/2). Also assume that the probability of presence in each of the finer domains is p1 and p2. Assuming independence between the two fine-grain domains, the probability of presence at the coarse domain will be 1 – (1 – p1)(1 – p2). For a simple case, assume that p1 = p2:

Probability of occurrence at coarse scale as a function of probability of occurrence at fine scales

We can see that except at the extreme cases of p1 = p2 = 0 and 1, coarse-scale probability of occurrence is always higher than fine-scale probability of occurrence. So the relevant question is “Does the increased probability of occurrence at coarse scales exceed what we’d expect by chance given that the coarse domain is composed of fine domains?” If so, only then can we say that temperature is a more important determinant of distribution at coarse scales. And that is what I would take as verification of this particular prediction of the Eltonian noise hypothesis.

Citation

Schweiger, A.H. and Beierkuhnlein, C. 2016. Scale dependence of temperature as an abiotic driver of species’ distributions. Global Ecology and Biogeography 25:1013-1021. DOI: 10.1111/geb.12463

Posted in biogeographic range, grain size, Nichey, spatial scale, temperature | No Comments »

Species distribution models not for species

July 20th, 2016

Mathematically these are all the same.

Have you ever see the number of people who drowned by falling into a swimming pool–films starring Nicholas Cage model? You might also know it as “linear regression.” Have you ever seen a species distribution model? By calling it thus we make the same limiting semantic complexification as in the first case.

This post is not about the debate over whether we should be calling it species distribution modeling or ecological niche modeling (notice the participle form of each term–adding an ing refers to the act of modeling). I’m talking about calling them species distribution models.

To put it bluntly, the underlying mathematics of a SDM doesn’t care whether it’s depicting a species or anything else. In fact there are numerous examples of applications of species-less “SDMs”:

Snakebite prevalence (Yañez-Arenas et al. 2016).
Presence of Pygmies in Africa (Olivero et al. 2016)
Safe sites for rescue helicopters to land in Yosemite National Park (Doherty et al. 2013)

I even once met a person who uses Maxent to locate opportune places for underwater archaeology!

The fundamental commonality that allows all of these phenomena to be modeled by “SDM” algorithms is the nature of the response data–either unary (i.e., just presences) or binary (presences and absences). Ergo, if you can describe a pattern with unary or binary data, you can also likely apply an “SDM”.

So feel free to refer to “my species distribution model” to reference your particular model of a species’ distribution. But don’t let the moniker box you into thinking they’re just for species!

Posted in Nichey, species distribution modeling | No Comments »

The fine limits of climate data

July 13th, 2016

Modeling Geum radiatum was very challenging compared to modeling I’ve done before because we know the species specializes on small habitats with microclimates. And these particular microhabitats are are not well reflected by the coarse-scale climate data that is available. I used ClimateNA for the basic climate layers. ClimateNA purports to be “scale free” but actually isn’t (and its creators acknowledge this). Indeed, I found that the relative humidity measured by Eric’s team was far higher than the values predicted by ClimateNA for the same years and locations. Ergo, ClimateNA (and other similar climate data sets) cannot really reflect microclimate on the order of of 100s of meters or less.

So how can one capture microclimate when one doesn’t have microclimate rasters? One option is to use microtopography, fine-scale measures of slope and aspect that hopefully help create microclimatic conditions relevant to the species. This gives us the “topoclimate” of Guitiérrez Illán et al. (2010). So, I coupled fine-scale topography data at 90-meter resolution with climate data from ClimateNA at 1-km resolution.

Now the normal way to combine data sets of differing resolutions would be to resample the coarser-scaled rasters so they have the same resolution as the fine-scaled rasters. But I felt this would create an illusion of precision in how well we know climate at fine scales.

So instead I wrote an R script to extract values from each data set in its native resolution and fed that to the predict() function. Essentially, the script goes cell-by-cell across the fine-scaled topography raster stack and extracts the matching coarse-scaled data from the climate raster stack. Then it writes a fine-scaled prediction raster from the ENM.

Geum radiatum (zoomed prediction)

The output doesn’t look especially blocky when you see it in the publication (Ulrey et al. 2016), but if you look carefully at the zoom-ins you can indeed see the effect of matching coarse- and fine-scaled data sets. The larger reddish cells reflect appropriate macroclimatic conditions (i.e., at 1000-m scales), while the smaller dark red areas indicate appropriate macro- and topoclimatic conditions.

I admit the map isn’t visually pleasing–it’s obviously an artifact that the predictions are blocky. But by keeping the data sets at their native resolution I believe the output better reflects the limit of our knowledge about the predictors, and in turn, the geographic realization of Geum‘s microclimatic niche.

Citations

Guitiérrez Illán, J., Gutiérrez, D., and Wilson, R.J. 2010. The contributions of topoclimate and land cover to species distributions and abundance: Fine-resolution tests for a mountain butterfly fauna. Global Ecology and Biogeography 19:159-173. (article page)

Ulrey, C., Quiantana-Ascencio, P.F., Kauffman, G., Smith, A.B., and Menges, E.S. 2016. Life at the top: Long-term demography, microclimatic refugia, and responses to climate change for a high-elevation southern Appalachian endemic plant. Biological Conservation 200:80-92. (article page)

Tags: climate, macroclimate, microclimate, spatial scale, topoclimate
Posted in Nichey, species distribution modeling | No Comments »

Cliffhanger

July 8th, 2016

I think I was on a long-haul flight across the Pacific when I succumbed to jet-lag induced doldrums and watched Sylvester Stallone’s Cliffhanger which stars him (surprise) as a mountaineer who gets himself out of a dastardly plot by climbing around and flexing his muscles. So if there’s a Rocky of the rare plant world, it’s Appalachian avens, or Geum radiatum.

Like Stallone, Geum likes to hang on cliffs in the Southern Appalachians. And talk about hang! Eric Menges sent me a few photos of their field work–they use ladders to census the populations. Now there are also a few populations found on so-called grassy balds (I am fond of this name), but the cliffside populations tend to be more common. And these cliffs occur at high elevations where the distinction between cloud and mist dissipates.

Geum radiatum (Wikimedia)

And like Stallone Geum is threatened–so much so it’s listed on the US Endangered Species Act. In fact Pedro Quintana-Ascencio‘s models predict that the overall growth rate of the populations for which they have census data is currently below replacement level. But–intriguingly–growth rate is highly positively linked to relative humidity.

Which is really cool in a technical way–because I was able to extract coarse-scale humidity data from the ClimateNA dataset, then relate this to their fine-scale measurements. And this in turn allowed Pedro to predict population growth rates under future scenarios of climate change.

“Hang in there, Geum!”

Alas, Geum needs a strongman like Stallone… even though we predicted relative humidity will drop by just a few percent, that’s enough to exacerbate the species’ current downward trajectory. It’s a slow fall, but it’s still a fall.

Citations

Tags: climate change, climate change vulnerabilty, endangered species, Geum radiatum, microclimate, rare species
Posted in global change, Nichey, species distribution modeling | No Comments »