OHI Science | Citation policy
[REFERENCE RMD FILE: https://cdn.rawgit.com/OHI-Science/ohiprep/master/globalprep/hab_seaice/v2017/hab_seaice_dataprep.html]
This calls a series of scripts to download and analyze sea ice data for the habitat subgoal and coastal protection goal.
An additional year of data has been added. The raster used to extract the region data is now created using the fasterize package which we have found to be more reliable. We attempted to transition to using the sf package for the spatial data, however, this was mostly not possible due to our reliance on the raster::extract function which is not supported for the sf format. Maybe next year we can make this update.
#Data Source Reference:
Cavalieri, D.J., Parkinson, C.L., Gloersen, P. and Zwally, H. (2014). Sea ice concentrations from Nimbus-7 SMMR and DMSP SMM/I-SSMIS passive microwave data. 1979-2014. Boulder, Colorado, USA: NASA National Snow and Ice Data Center Distributed Active Archive Center (https://doi.org/10.5067/8GQ8LZQVL0VL).
Downloaded: September 19 2017
Description: Monthly sea ice extent data.
Data are in the following FTP directory: ftp://sidads.colorado.edu/pub/DATASETS/nsidc0051_gsfc_nasateam_seaice/final-gsfc/
Within the final-gsfc directory are north and south directories that contain data files, and a browse directory that contains browse image files. Daily and monthly data are further separated into directories named daily and monthly. For final daily data, there is also one directory for each year of available data. For example, all of the north daily data for 1990 are in a directory named /nsidc0051_gsfc_nasateam_seaice/final-gsfc/north/daily/1990/.
The directory structure is illustrated below; not all directory levels are shown. The structure for each south directory matches that of the corresponding north directory. Each browse directory is divided into a structure that reflects that of the data. In this illustration, the year directories underneath final-gsfc are representative placeholders; there are actually many such directories, each named for the year of data it contains, such as 1987, 2000, etc.
/nsidc0051_gsfc_nasateam_seaice . . /final-gsfc . . . . /browse . . . . . . /north . . . . . . . . /daily . . . . . . . . . . /year . . . . . . . . /monthly . . . . . . /south . . . . /north . . . . . . /daily . . . . . . . . /year . . . . . . /monthly . . . . /south
For complete documentation and more information about data access, please see:
http://nsidc.org/data/nsidc-0051.html
If you wish to be notified of updates or corrections to these data, please register with NSIDC User Services by sending e-mail to:
nsidc@nsidc.orgIdentify yourself as a user of “Sea Ice Concentrations from Nimbus-7 SMMR and DMSP SSM/I-SSMIS Passive Microwave Data (NSIDC-0051).” Include your name, e-mail address, postal address, and telephone number.
If you have questions, please contact NSIDC User Services.
CONTACT INFORMATION: User Services National Snow and Ice Data Center CIRES, 449 UCB University of Colorado Boulder, CO USA 80309-0449 Phone: +1 303-492-6199 Fax: +1 303-492-2468 E-mail: nsidc@nsidc.org
Time range: 1970-2015
# load libraries, set directories
library(sp)
library(raster)
library(fasterize)##
## Attaching package: 'fasterize'## The following object is masked from 'package:graphics':
##
## plotlibrary(sf)## Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.8.0, lwgeom 2.2.2 r14797library(rgdal)## rgdal: version: 1.2-15, (SVN revision 691)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 2.1.0, released 2016/04/25
## Path to GDAL shared files: /usr/share/gdal/2.1
## GDAL binary built with GEOS: TRUE
## Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]
## Path to PROJ.4 shared files: (autodetected)
## Linking to sp version: 1.2-5library(fields) #colors in Status_Trend.R## Loading required package: spam## Loading required package: dotCall64## Loading required package: grid## Spam version 2.1-1 (2017-07-02) is loaded.
## Type 'help( Spam)' or 'demo( spam)' for a short introduction
## and overview of this package.
## Help for individual functions is also obtained by adding the
## suffix '.spam' to the function name, e.g. 'help( chol.spam)'.##
## Attaching package: 'spam'## The following objects are masked from 'package:base':
##
## backsolve, forwardsolve## Loading required package: maps## comment out when knitting
#setwd("globalprep/hab_seaice/v2017")
source('../../../src/R/common.R') # directory locationsThese maps were made by the PreparingSpatialFiles.R script. If there are changes to the OHI regions, this script needs to be run. Additionally, it is critical to walk through the ObtainingData.R script if any of the spatial files have been modified (this saves the files as spatial points).
The original polygon files are used from: git-annex/globalprep/_raw_data/NSIDC_SeaIce/v2015 However, we use the new fasterize package to rasterize the polygons and save them to: git-annex/globalprep/_raw_data/NSIDC_SeaIce/v2017
maps <- file.path(dir_M, "git-annex/globalprep/_raw_data/NSIDC_SeaIce/v2017")# identify the year to save raw data file (Need to create this file)
assessYear <- 'v2017'
# final year of data:
final.year <- 2016
# Establish: CRS, website to collect data, data selection parameters
pixel = 25000 # pixel dimension in meters for both x and y
# epsg projection 3411 - nsidc sea ice polar stereographic north (http://spatialreference.org/ref/epsg/3411/)
# epsg projection 3412 - nsidc sea ice polar stereographic south (http://spatialreference.org/ref/epsg/3412/)
prj.n = '+proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs'
prj.s = '+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs'
prj.mol = '+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs'
# URL base (ub), filename format for final monthly data is nt_YYYYMM_SSS_vVV_R.bin
ub.n = 'ftp://sidads.colorado.edu/pub/DATASETS/nsidc0051_gsfc_nasateam_seaice/final-gsfc/north/monthly'
ub.s = 'ftp://sidads.colorado.edu/pub/DATASETS/nsidc0051_gsfc_nasateam_seaice/final-gsfc/south/monthly'
poles = c('n','s')
years = c(1979:final.year) #Full range of data
months = 1:12
n.pym=length(poles)*length(years)*length(months)
i.pym = 0
t0 = Sys.time()Collects the data for each month/year from the website and add to raster stack that is saved in tmp folder as: n_rasters_points.rdata or s_rasters_points.rdata. And, if it doesn’t already exist, it converts the region shapefile into a raster points file.
source("ObtainingData.R")Using the data from the .rdata files created in Function 1, this function calculates status and trend for shoreline ice and ice edge habitat. Data is saved in intermediate folder:
ref.years <- 1979:2000
source("Status_Trend.R")n_edge <- read.csv("int/n_IceEdgeHabitat_ref1979to2000.csv")
s_edge <- read.csv("int/s_IceEdgeHabitat_ref1979to2000.csv")
edge <- rbind(n_edge, s_edge)
edge <- edge %>%
dplyr::filter(Reference_avg1979to2000monthlypixels != 0) %>%
dplyr::filter(!(rgn_id %in% c(59, 141, 219, 4, 172, 94))) %>% #anomolous eez regions with very little ice cover
dplyr::filter(!(rgn_id %in% c("248300", "258510", "258520", "258600", "258700"))) %>% # ccamlr: cut some regions due to minimal ice (<200 km2 per year - average of months)
dplyr::filter(rgn_nam != "DISPUTED") %>%
dplyr::mutate(habitat="seaice_edge")
n_shore <- read.csv("int/n_IceShoreProtection_ref1979to2000.csv")
s_shore <- read.csv("int/s_IceShoreProtection_ref1979to2000.csv")
shore <- rbind(n_shore, s_shore)
shore <- shore %>%
dplyr::filter(Reference_avg1979to2000monthlypixels != 0) %>%
dplyr::filter(!(rgn_id %in% c(59, 89, 177, 178))) %>% #anomolous eez regions with very little ice cover
dplyr::filter(rgn_nam != "DISPUTED") %>%
dplyr::mutate(habitat="seaice_shoreline")
data <- rbind(edge, shore)
data <- data %>%
dplyr::mutate(km2 = Reference_avg1979to2000monthlypixels/12 * (pixel/1000)^2)
write.csv(data, "int/sea_ice.csv", row.names=FALSE)
## health data
health <- data %>%
dplyr::filter(rgn_typ == "eez") %>%
dplyr::select(rgn_id, habitat, dplyr::starts_with('pctdevR')) %>%
tidyr::gather("year", "health", -(1:2)) %>%
dplyr::mutate(year = substring(year, 9, 12)) %>%
dplyr::mutate(year = as.numeric(year)) %>%
dplyr::mutate(health = ifelse(health > 1, 1, health))
write.csv(health, "output/hab_ice_health_eez.csv", row.names=FALSE)
## trend data
trend <- data %>%
dplyr::filter(rgn_typ == "eez") %>%
dplyr::select(rgn_id, habitat, dplyr::starts_with('Trend')) %>%
tidyr::gather("year", "trend", -(1:2)) %>%
dplyr::mutate(year = substring(year, 13, 16)) %>%
dplyr::mutate(year = as.numeric(year)) %>%
dplyr::mutate(trend = trend * 5) %>%
dplyr::mutate(trend = ifelse(trend>1, 1, trend)) %>%
dplyr::mutate(trend = ifelse(trend<(-1), -1, trend))
write.csv(trend, "output/hab_ice_trend_eez.csv", row.names=FALSE)
### get extent data
extent <- data %>%
dplyr::filter(rgn_typ == "eez") %>%
dplyr::mutate(year = 2016) %>%
dplyr::select(rgn_id, habitat, year, km2)
write.csv(extent, "output/hab_ice_extent_eez.csv", row.names=FALSE)Compare to last year’s data. There should be a strong correlation between these data. This is what we observed in the data.
Pair with region names for sanity check.
## Health
health <- read.csv("output/hab_ice_health_eez.csv")
# compare to last year's data to make sure everything looks basically fine
old <- read.csv("../v2016/output/hab_ice_health_eez_updated.csv") %>%
dplyr::rename(old_health = health) %>%
dplyr::left_join(health, by = c("rgn_id", "habitat", "year"))
plot(old$old_health, old$health)
abline(0, 1, col="red")
## Trend
trend <- read.csv("output/hab_ice_trend_eez.csv")
old <- read.csv("../v2016/output/hab_ice_trend_eez_updated.csv") %>%
dplyr::rename(old_trend = trend) %>%
dplyr::left_join(trend, by = c("rgn_id", "habitat", "year"))
plot(old$old_trend, old$trend)
abline(0, 1, col="red")
## Extent
# compare to last year's data to make sure everything looks basically fine
extent <- read.csv("output/hab_ice_extent_eez.csv")
old <- read.csv("../v2016/output/hab_ice_extent_eez_updated.csv") %>%
dplyr::rename(old_km2 = km2) %>%
dplyr::left_join(extent, by = c("rgn_id", "year", "habitat"))
plot(old$old_km2, old$km2)
abline(0, 1, col="red")
## Pair with region names
library(ohicore)
regions <- rgn_master %>%
dplyr::select(rgn_id = rgn_id_2013, rgn_name = rgn_nam_2013) %>%
unique()
data <- read.csv('output/hab_ice_health_eez.csv') %>%
dplyr::left_join(regions, by="rgn_id") %>%
dplyr::arrange(habitat, health)