Base run

A simple run of xagg, aggregating gridded temperature data over US counties. For a deeper dive into xagg’s functionality, see the Detailed Code Run.

import xagg as xa
import xarray as xr
import numpy as np
import geopandas as gpd

Import

The sample data in this example are:

• gridded: month-of-year average temperature projections for the end-of-century from a climate model (CCSM4)

• shapefiles: US counties

# Load some climate data as an xarray dataset
ds = xr.open_dataset('../../data/climate_data/tas_Amon_CCSM4_rcp85_monthavg_20700101-20991231.nc')
ds

<xarray.Dataset>
Dimensions:   (lon: 288, lat: 192, month: 12, bnds: 2)
Coordinates:
    height    float64 ...
  * lon       (lon) float64 0.0 1.25 2.5 3.75 5.0 ... 355.0 356.2 357.5 358.8
  * lat       (lat) float64 -90.0 -89.06 -88.12 -87.17 ... 88.12 89.06 90.0
  * month     (month) int64 1 2 3 4 5 6 7 8 9 10 11 12
Dimensions without coordinates: bnds
Data variables:
    lat_bnds  (month, lat, bnds) float64 ...
    lon_bnds  (month, lon, bnds) float64 ...
    tas       (month, lat, lon) float32 ...

xarray.Dataset

• Dimensions:
  • lon: 288
  • lat: 192
  • month: 12
  • bnds: 2
• Coordinates: (4)
  • height
    ()
    float64
    ...

    units :

    m

    axis :

    Z

    positive :

    up

    long_name :

    height

    standard_name :

    height

    [1 values with dtype=float64]

  • lon
    (lon)
    float64
    0.0 1.25 2.5 ... 356.2 357.5 358.8

    bounds :

    lon_bnds

    units :

    degrees_east

    axis :

    X

    long_name :

    longitude

    standard_name :

    longitude

    array([  0.  ,   1.25,   2.5 , ..., 356.25, 357.5 , 358.75])

  • lat
    (lat)
    float64
    -90.0 -89.06 -88.12 ... 89.06 90.0

    bounds :

    lat_bnds

    units :

    degrees_north

    axis :

    Y

    long_name :

    latitude

    standard_name :

    latitude

    array([-90.      , -89.057592, -88.115183, -87.172775, -86.230366, -85.287958,
           -84.34555 , -83.403141, -82.460733, -81.518325, -80.575916, -79.633508,
           -78.691099, -77.748691, -76.806283, -75.863874, -74.921466, -73.979058,
           -73.036649, -72.094241, -71.151832, -70.209424, -69.267016, -68.324607,
           -67.382199, -66.439791, -65.497382, -64.554974, -63.612565, -62.670157,
           -61.727749, -60.78534 , -59.842932, -58.900524, -57.958115, -57.015707,
           -56.073298, -55.13089 , -54.188482, -53.246073, -52.303665, -51.361257,
           -50.418848, -49.47644 , -48.534031, -47.591623, -46.649215, -45.706806,
           -44.764398, -43.82199 , -42.879581, -41.937173, -40.994764, -40.052356,
           -39.109948, -38.167539, -37.225131, -36.282723, -35.340314, -34.397906,
           -33.455497, -32.513089, -31.570681, -30.628272, -29.685864, -28.743455,
           -27.801047, -26.858639, -25.91623 , -24.973822, -24.031414, -23.089005,
           -22.146597, -21.204188, -20.26178 , -19.319372, -18.376963, -17.434555,
           -16.492147, -15.549738, -14.60733 , -13.664921, -12.722513, -11.780105,
           -10.837696,  -9.895288,  -8.95288 ,  -8.010471,  -7.068063,  -6.125654,
            -5.183246,  -4.240838,  -3.298429,  -2.356021,  -1.413613,  -0.471204,
             0.471204,   1.413613,   2.356021,   3.298429,   4.240838,   5.183246,
             6.125654,   7.068063,   8.010471,   8.95288 ,   9.895288,  10.837696,
            11.780105,  12.722513,  13.664921,  14.60733 ,  15.549738,  16.492147,
            17.434555,  18.376963,  19.319372,  20.26178 ,  21.204188,  22.146597,
            23.089005,  24.031414,  24.973822,  25.91623 ,  26.858639,  27.801047,
            28.743455,  29.685864,  30.628272,  31.570681,  32.513089,  33.455497,
            34.397906,  35.340314,  36.282723,  37.225131,  38.167539,  39.109948,
            40.052356,  40.994764,  41.937173,  42.879581,  43.82199 ,  44.764398,
            45.706806,  46.649215,  47.591623,  48.534031,  49.47644 ,  50.418848,
            51.361257,  52.303665,  53.246073,  54.188482,  55.13089 ,  56.073298,
            57.015707,  57.958115,  58.900524,  59.842932,  60.78534 ,  61.727749,
            62.670157,  63.612565,  64.554974,  65.497382,  66.439791,  67.382199,
            68.324607,  69.267016,  70.209424,  71.151832,  72.094241,  73.036649,
            73.979058,  74.921466,  75.863874,  76.806283,  77.748691,  78.691099,
            79.633508,  80.575916,  81.518325,  82.460733,  83.403141,  84.34555 ,
            85.287958,  86.230366,  87.172775,  88.115183,  89.057592,  90.      ])

  • month
    (month)
    int64
    1 2 3 4 5 6 7 8 9 10 11 12

    array([ 1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12])

• Data variables: (3)
  • lat_bnds
    (month, lat, bnds)
    float64
    ...

    [4608 values with dtype=float64]

  • lon_bnds
    (month, lon, bnds)
    float64
    ...

    [6912 values with dtype=float64]

  • tas
    (month, lat, lon)
    float32
    ...

    [663552 values with dtype=float32]

• Indexes: (3)
  • lon
    PandasIndex

    PandasIndex(Index([   0.0,   1.25,    2.5,   3.75,    5.0,   6.25,    7.5,   8.75,   10.0,
            11.25,
           ...
            347.5, 348.75,  350.0, 351.25,  352.5, 353.75,  355.0, 356.25,  357.5,
           358.75],
          dtype='float64', name='lon', length=288))

  • lat
    PandasIndex

    PandasIndex(Index([            -90.0, -89.0575916230366, -88.1151832460733,
           -87.1727748691099, -86.2303664921466, -85.2879581151832,
           -84.3455497382199, -83.4031413612565, -82.4607329842932,
           -81.5183246073298,
           ...
            81.5183246073298,  82.4607329842932,  83.4031413612565,
            84.3455497382199,  85.2879581151832,  86.2303664921466,
            87.1727748691099,  88.1151832460733,  89.0575916230366,
                        90.0],
          dtype='float64', name='lat', length=192))

  • month
    PandasIndex

    PandasIndex(Index([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12], dtype='int64', name='month'))

• Attributes: (0)

# Load US counties shapefile as a geopandas GeoDataFrame
gdf = gpd.read_file('../../data/geo_data/UScounties.shp')
gdf.head()

	NAME	STATE_NAME	STATE_FIPS	CNTY_FIPS	FIPS	geometry
0	Lake of the Woods	Minnesota	27	077	27077	POLYGON ((-95.34283 48.54668, -95.34105 48.715...
1	Ferry	Washington	53	019	53019	POLYGON ((-118.85163 47.94956, -118.84846 48.4...
2	Stevens	Washington	53	065	53065	POLYGON ((-117.43883 48.04412, -117.54219 48.0...
3	Okanogan	Washington	53	047	53047	POLYGON ((-118.97209 47.93915, -118.97406 47.9...
4	Pend Oreille	Washington	53	051	53051	POLYGON ((-117.43858 48.99992, -117.03205 48.9...

Aggregate

Now, aggregate the gridded variable in ds onto the polygons in gdf. Use the option silent=True if you’d like to suppress the status updates, or set it as a default using xa.set_defaults(silent=True).

# Calculate overlaps
weightmap = xa.pixel_overlaps(ds,gdf)

creating polygons for each pixel...
calculating overlaps between pixels and output polygons...
success!

# Aggregate
aggregated = xa.aggregate(ds,weightmap)

adjusting grid... (this may happen because only a subset of pixels were used for aggregation for efficiency - i.e. [subset_bbox=True] in xa.pixel_overlaps())
grid adjustment successful
aggregating tas...
all variables aggregated to polygons!

Convert

Finally, convert the aggregated data back into the format you would like.

# Example as an xarray dataset
ds_out = aggregated.to_dataset()
ds_out

<xarray.Dataset>
Dimensions:     (poly_idx: 3141, month: 12)
Coordinates:
  * poly_idx    (poly_idx) int64 0 1 2 3 4 5 6 ... 3135 3136 3137 3138 3139 3140
  * month       (month) int64 1 2 3 4 5 6 7 8 9 10 11 12
Data variables:
    NAME        (poly_idx) object 'Lake of the Woods' 'Ferry' ... 'Broomfield'
    STATE_NAME  (poly_idx) object 'Minnesota' 'Washington' ... 'Colorado'
    STATE_FIPS  (poly_idx) object '27' '53' '53' '53' ... '02' '02' '02' '08'
    CNTY_FIPS   (poly_idx) object '077' '019' '065' '047' ... '240' '068' '014'
    FIPS        (poly_idx) object '27077' '53019' '53065' ... '02068' '08014'
    tas         (poly_idx, month) float64 263.9 268.8 274.0 ... 276.4 270.4

xarray.Dataset

• Dimensions:
  • poly_idx: 3141
  • month: 12
• Coordinates: (2)
  • poly_idx
    (poly_idx)
    int64
    0 1 2 3 4 ... 3137 3138 3139 3140

    array([   0,    1,    2, ..., 3138, 3139, 3140])

  • month
    (month)
    int64
    1 2 3 4 5 6 7 8 9 10 11 12

    array([ 1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12])

• Data variables: (6)
  • NAME
    (poly_idx)
    object
    'Lake of the Woods' ... 'Broomfi...

    array(['Lake of the Woods', 'Ferry', 'Stevens', ...,
           'Southeast Fairbanks', 'Denali', 'Broomfield'], dtype=object)

  • STATE_NAME
    (poly_idx)
    object
    'Minnesota' ... 'Colorado'

    array(['Minnesota', 'Washington', 'Washington', ..., 'Alaska', 'Alaska',
           'Colorado'], dtype=object)

  • STATE_FIPS
    (poly_idx)
    object
    '27' '53' '53' ... '02' '02' '08'

    array(['27', '53', '53', ..., '02', '02', '08'], dtype=object)

  • CNTY_FIPS
    (poly_idx)
    object
    '077' '019' '065' ... '068' '014'

    array(['077', '019', '065', ..., '240', '068', '014'], dtype=object)

  • FIPS
    (poly_idx)
    object
    '27077' '53019' ... '02068' '08014'

    array(['27077', '53019', '53065', ..., '02240', '02068', '08014'],
          dtype=object)

  • tas
    (poly_idx, month)
    float64
    263.9 268.8 274.0 ... 276.4 270.4

    array([[263.91894338, 268.83407312, 273.97753272, ..., 283.79866008,
            275.10910005, 266.01617602],
           [271.79416912, 275.6313644 , 276.94708001, ..., 281.64892744,
            276.72788621, 272.2569337 ],
           [272.11315505, 275.91027861, 277.35535375, ..., 282.07634405,
            277.01922179, 272.51605575],
           ...,
           [262.8463124 , 263.0001854 , 265.4380367 , ..., 271.71306123,
            264.94652565, 263.29793632],
           [265.08434154, 264.54793604, 267.20395358, ..., 272.51411657,
            266.14508751, 265.68265969],
           [270.80386353, 273.4302063 , 275.95550537, ..., 283.54470825,
            276.38360596, 270.44485474]])

• Indexes: (2)
  • poly_idx
    PandasIndex

    PandasIndex(Index([   0,    1,    2,    3,    4,    5,    6,    7,    8,    9,
           ...
           3131, 3132, 3133, 3134, 3135, 3136, 3137, 3138, 3139, 3140],
          dtype='int64', name='poly_idx', length=3141))

  • month
    PandasIndex

    PandasIndex(Index([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12], dtype='int64', name='month'))

• Attributes: (0)

# Example as a pandas dataframe
df_out = aggregated.to_dataframe()
df_out

		NAME	STATE_NAME	STATE_FIPS	CNTY_FIPS	FIPS	tas
poly_idx	month
0	1	Lake of the Woods	Minnesota	27	077	27077	263.918943
	2	Lake of the Woods	Minnesota	27	077	27077	268.834073
	3	Lake of the Woods	Minnesota	27	077	27077	273.977533
	4	Lake of the Woods	Minnesota	27	077	27077	283.141960
	5	Lake of the Woods	Minnesota	27	077	27077	290.623952
...	...	...	...	...	...	...	...
3140	8	Broomfield	Colorado	08	014	08014	297.646820
	9	Broomfield	Colorado	08	014	08014	292.368988
	10	Broomfield	Colorado	08	014	08014	283.544708
	11	Broomfield	Colorado	08	014	08014	276.383606
	12	Broomfield	Colorado	08	014	08014	270.444855

37692 rows × 6 columns