Climate Change in 4 Dimensions

From April 8 2104 till June 17 I participated in online course

www.coursera.org

Climate Change in 4 Dimensions
University of San Diego.

Overal I learnt a lot about climate change beyond the science aspects

During the course I kept notes to provide some feedback and ideas for improvement.
I captured this in this post.

Lecture 5

Early part is very US centric 

Global Average Warming Chart after Ocean Acidification and Coral , the text is to small to read.

Lecture 6

Very informative and enjoyable, Professor David G. Victor engaged me well.Introduced Stock Problem Economic model to me. Already introducing new ideas to my thoughts

Lecture 9

Weekly activity
Survey is USA oriented. Of course the US should act independent of other countries. The question for me is should Australia act independent of other countries.

This weeks graded quiz was not up to scratch. Usually only one or none of the question are badly worded, easily to be misunderstood, or difficult to attribute to any of the reading or lectures. This week multiple questions including.
Question 8 , 9. difficult to understand/read. 
"Heat waves and heavy precipitation" is not a trend.    

Lecture 10

Feedback: 
Slide from Rogell, Meinshausen 2012 is out of date, wrong name inconsistent with other slides. 
Chaparral is not a common international word, and is not what would grow in Australia. 
Question in the ungraded quiz unclear

Based on global climate models, what regions will see the biggest precipitation increase with global warming?

northern

western

southern

eastern

This question is about North America, Not clear

Concepts

Uses indices to show warming particular in USA, definitely regional/spatial
Indices

• Mean annual temperature
• Number of cold days
• Number of hot days
• Agricultural Based - frost, super hot days, frost free (growing season) forecast.
• Keeling curve shows growing season has increased by about 10 days, in agreement with ag-models

Heat Wave

• Multiple factors some of which more common due to CC
• Europe 2003, 2006, Australia 2009 (Karoly), Russia 2010, USA 2009

Models don't seem to capture very cold days. 

Midterm

Question 7 unclear
If someone makes a prediction, and it comes true, does that mean that the prediction is correct?
does it mean hypothesis or prediction is true?

Question 9 
Debatable if one or more answers are correct

Week 8 

Lecture 14 
Ice, snow and water
Does NOT mention ice West antarctica Ice Shelf info?
Projections continue to change with increasing understanding?
Parts of the lecture are dating quickly.

Lecture 15
Very directed to a Californian Audience. The level of knowledge about geography of California was assumed to be higher than my own, though I have travelled there.  Think this would exclude some of the coursera audience.

Quiz  Question 9

"Please fill in the blank: By 2050, the number of ‘extremely hot’ days could increase ________."

This statement was not in the lecture notes, so how does one know what location of study to cite? Does it it relate to global, USA, California, San Diego, or Australia extremely hot days.?

Lecture 19

Spelling mistake in the "Check your knowledge quiz" question 4
Arpanet

Thoughts on Pattern Scaling

During the last week I have been consider ideas put forward at the

Pattern Scaling, Climate Model Emulators and their Application to the new Scenario Process
NCAR, Boulder Colorado, April 23-25 2014

and

Lopez et. la. 2013 Robustness of pattern scaled climate change scenarios for adaption decision support

Major Aims

• "Fit empirical / statistical relation b/w impact relevant climate variables and large scale quantities obtainable through simple models"
• "Run simple models under arbitrary scenarios and recover impact relevant outcomes by applying those relations"

• User Needs
• Impact research
• policy makers
• Social Economic
• higher resolution
• 
  
• Uncertainty
• Handling
• Quantitating 
• sufficiently low uncertainty for outcome information produced to be useful.

Standard Pattern Scaling

• Developed, tested and applied for 20 years
• provide a simplified representation of climate system responses.
• 
  
  • local (or regional) changes in these variables tend to increase linearly with the global warming over the coming century.
  1. local change can be seen as a 'response' to the global warming (GW)
  
  critical assumption is the there is linear relationship b/w a scalar, and a geographical response pattern

Flaws / Concerns

• main climate mechanisms are not linear
• feedback
• timescale in response change
• patterns evolve
• 
  

Uncertainty 

Uncertainty hard to capture

• model uncertainty
• Multimodel ensemble can reduce uncertainty
• scenario uncertainty
• depend on statistical assumption
• analysis of variance 
• map std dev.

NCL Regional Temperature Anomolies

Data
HadCRUT 4
Surface Temperature Anomalies (C with respect to 1961-1990)

Choose a region

Getting my regions in a reusable, correctly projected, and regridded  was more difficult than I thought it would be.
Publicly available global region based netCDF files are not common
In the end, I used a shape file to create a netCDF with a single time dimension, and  integer layer of region ids at each 1 degree

Australia

Just going with continents at this stage and Australia is the country I was born.

Start with a single day.
; read only desired  time 1970
 x  = fin->temperature_anomaly(1440,:,:)

 xr = where(region.eq.2,x,x@_FillValue)

 print(avg(xr))

Then time series for a single year. Oct 2012 to Sept 2013 

; Time Series average for Australia 2013
 x1=fin->temperature_anomaly(1953:1964,:,:)
 rconform = conform(x1, region, (/1,2/))
 xr1 =mask(x1, rconform, 2)

 xa1 = dim_avg_n(xr1,(/1 , 2/))
 t1 = ispan(0,12,1)

 wks   = gsn_open_wks ("x11","xy")                ; open workstation
 res                  = True                     ; plot mods desired
 res@tiMainString     = "Average Australian Temperature Anomaly 2013" ;

Result similar to BOM 

Then annual time series for a decade or so

 yStart = 1990
 yEnd = 2013
 tStart  = (yStart - T_OFFSET ) * 12;
 tEnd = (yEnd - T_OFFSET) * 12 + 11
 x2=fin->temperature_anomaly(tStart:tEnd,:,:)
    
 rconform2 = conform(x2, region, (/1,2/))
 xr2 =mask(x2, rconform2, 2)
 copy_VarCoords(x2, xr2) ; need dim metadata retained for clim functions
 xa2 = dim_avg_n(xr2,(/1 , 2/))
  xannual = month_to_annual(xa2, 1)  ;Annual Average Temperautre             

  printVarSummary(xannual)  

 wks   = gsn_open_wks ("ps","xy")                ; open workstation
 res                  = True                     ; plot mods desired
 res@tiMainString     = "Annual Mean  Australian Temperature  1990-2013" ;
 res@tiYAxisString = "Anomalies" ; y-axis label      

 res@gsnYRefLine           = 0.              ; reference line   
 res@gsnXYBarChart         = True            ; create bar chart 
 res@gsnAboveYRefLineColor = "red"           ; above ref line fill red
 res@gsnBelowYRefLineColor = "blue"          ; below ref line fill blue

res@tiXAxisString = "Year"

 plot  = gsn_csm_xy (wks,ispan(yStart,yEnd,1),xannual,res) ; create plot  

Then summer time series for a decade

 ; Summer Time Average Temperature anomaly for 2000's
 yStart = 2000
 yEnd = 2013
 tStart  = (yStart - T_OFFSET ) * 12 ;
 tEnd = (yEnd - T_OFFSET) * 12 + 11
 x3=fin->temperature_anomaly(tStart:tEnd,:,:)
    
 rconform3 = conform(x3, region, (/1,2/)) ;Create Mask grid with time dim of data
 xr3 =mask(x3, rconform3, 2) ; mask out data not in region 2 , Australia
 copy_VarCoords(x3, xr3) ; need dim metadata retained for clim functions

 xseasonal = month_to_seasonN(xr3, (/ "DJF"/))  ; Seasonal Average Temperautre       
 printVarSummary(xseasonal)  

 xa3 = dim_avg_n(xseasonal,(/2 , 3/)) ; Average across  long and lat dimension
                                ; average across region as all other values masked 
 copy_VarCoords_2(xseasonal,xa3)
 printVarSummary(xa3)

 print(xa3)

 wks   = gsn_open_wks ("x11","xy")                ; open workstation
 res                  = True                     ; plot mods desired
 res@tiMainString     = "Summer Mean  Australian Temperature  1990-2013" ;
 res@tiYAxisString = "Anomalies" ; y-axis label      

 res@gsnYRefLine           = 0.              ; reference line   
 res@gsnXYBarChart         = True            ; create bar chart 
 res@gsnAboveYRefLineColor = "red"           ; above ref line fill red
 res@gsnBelowYRefLineColor = "blue"          ; below ref line fill blue

res@tiXAxisString = "Year"

 plot  = gsn_csm_xy (wks,ispan(yStart,yEnd,1),xa3(0,:),res) ; create plot 

Happy with this time to move on

Install CDO

Mac OSX Mavericks

Download https://code.zmaw.de/projects/cdo/files

In terminal

gunzip cdo-current.tar.gz

tar -xf cdo-current.tar

cd cdo-1.6.4r5

./configure

make

sudo make install

Projections and Gridding

Today I have been thinking about projections

CMIP5 is netCDF the projection is provided in the Metadata, resolution 0.5° x 0.5°

I want to perform polygon regions analysis.

I need the 'masking' of the CMIP5 to be efficient. Visualisation is not a concern.

What file format for storing polygons should I use ?

• projection 
• resolution

I have been looking at different programming languages, particularly R and NCL.

NCL is built to handle  netCDF data efficiently.

But I am having trouble grasping how I can use my polygons efficiently. The polygons are initially shapefiles. The example polygon code focuses on using polygons in visualisation. The example shapefile masking code seems very inefficient. 

Focusing on NCL

If I know CMIP5 resolution and projection I should be able to create netCDF file with a layer / a variable / or slot in array -  per region , with binary values.

Max and Min Latitude of the region could be used to reduce data extract from the netCDF

; read only desired area & times

x     = in->SST(tStrt:tLast,{latS:latN},{lonL:lonR})

Merge / AND / if / where

@FillValue is kind of like null for this variable, ignored by many functions.

  if(.not.isatt(data,"_FillValue")) then
    data@_FillValue = default_fillvalue(typeof(data))          ;-- make sure "data" has a missing value
  end if

x is CMIP5 netCDF data

regions is all my polygons

 x and regions should have the same dimensions

 xr = where(regions.eq.17, x, x@_FillValue)

Or once I have a single binary object per region

 x = in->SST(tStrt:tLast,{r17.latS:r17.latN},{r17.lonL:r17.lonR})

 xr = where(r17, x, x@_FillValue)

There is memory concern here. 

I am going with the concept by region. that is 

For each region

For each netCDF file

Open netCDF file
Populate variable (X) with netCDF data by region and time period
Close netCDF file
Perform any calculation which will optimise memory footprint 

x <- X

Delete X, keep derived data.


Still got the region shape file. Regrid to standard CMIP5 

Ahhh there is no standard grid ding CMIP5

Initially going for bi-linearly interpolation  to 1x1 rectilinear grid 

Why Bilinear

Focus initially land
Noaa using in this http://www.esrl.noaa.gov/psd/ipcc/ocn/details.html
Linear Interpolation had a good ranking in 

Mora, C., et al. (2013). "Biotic and human vulnerability to projected changes in ocean biogeochemistry over the 21st century." PLoS biology 11(10): e1001682.




Why Rectilinear grid



Simple
Too many people think in rectangles
I don't like how the areas are so different
Consider variation at later stage


Why 1 x 1


Because resolution should be reasonable with all the CMIP5 0.5 degree data

NCL Try outs

Installing NCL

Put this in my .bash_profile file:

# NCL additions

export NCARG_ROOT=/Volumes/Data/Users/fmacgill/Projects/ncl62/

export PATH=$NCARG_ROOT/bin:$PATH

Library Issue

Had issues with a library, solve by installing gcc47 through macports

http://www.ncl.ucar.edu/Download/macosx.shtml#libgomp

Put this in my .bash_profile file:

   export DYLD_FALLBACK_LIBRARY_PATH=/opt/local/lib/gcc47/

clmMonLLT

Decide to try out NCL with this example code.

https://www.ncl.ucar.edu/Applications/climo.shtml

Downloaded 

clim0_4.ncl

https://www.ncl.ucar.edu/Applications/Data/

Downloaded 

xieArkin_T42.nc

Want to see output immediately so changed Line 31

wks   = gsn_open_wks("x11" ,"climo")        ; open ps file

So what going on

Initial Region Masking map

There are many ways to break up the global into regions.
Already spent hours trying to find ways of classification of create regions
Already spent hours trying to locate shape files, geoTiffs, or any GIS compatible sources

Waste of time, the data exists, others have it. Need to work out who to ask and gain the confident to ask.

Need to focus on the process

Josh O'Brien
So from http://stackoverflow.com/questions/20146809/how-can-i-plot-a-continents-map-with-r
SRES are IPCC Emission Scenario Regions


Same code tweaked.

library(rworldmap)
library(rgeos)

sPDF <- getMap()
sres <-
    sapply(levels(sPDF$SRES),
           FUN = function(i) {
               ## Merge polygons within a continent
               poly <- gUnionCascaded(subset(sPDF, SRES==i))
               ## Give each polygon a unique ID
               poly <- spChFIDs(poly, i)
               ## Make SPDF from SpatialPolygons object
               SpatialPolygonsDataFrame(poly,
                                        data.frame(SRES=i, row.names=i))
           },
           USE.NAMES=TRUE)

## Bind the 11 SRES-level SPDFs into a single SPDF
sres <- Reduce(spRbind, sres)

## Plot to check that it worked
plot(sres, col=heat.colors(nrow(sres)))

## Check that it worked by looking at the SPDF's data.frame
## (to which you can add attributes you really want to plot on)
data.frame(sres)

Save for future use. Use  .RData for faster loading in the future.
> class(sres)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"
> proj4string(sres)
[1] "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"


> setwd("~/Projects/R/test")
> getwd()
[1] "/Volumes/Data/Users/fmacgill/Projects/R/test"

> sres <- data_name
> plot(sres, col=heat.colors(nrow(sres)))
> data_name <- sres
> save(data_name, file="SRESGlobal.RData" )

Clear Variables and check retrieval

> sres <- 1
> data_name <- 1
> load("SRESGlobal.RData")
> sres <- data_name
> plot(sres, col=heat.colors(nrow(sres)))