`scatstat1` documentation

The `scatstat1` function returns statistical values of all points within a given 1D radius of each value. This is similar to taking a moving mean, but points do not have to be equally spaced, nor do x values need to be monotonically increasing.

See also `scatstat2`.

Back to Climate Data Tools Contents.

## Contents

## Syntax

ybar = scatstat1(x,y,radius) ybar = scatstat1(x,y,radius,fun) ybar = scatstat1(...,options)

## Description

`ybar = scatstat1(x,y,radius)` returns the mean of all `y` values within specified `radius` at each point `x`.

`ybar = scatstat1(x,y,radius,fun)` applies any function `fun` to `y` values, default `fun` is `@mean`, but could be `@median`, `@nanstd`, etc.

`ybar = scatstat1(...,options)` allows any options that the function accepts. For example, `'omitnan'`.

## Example 1

Get the local median of all points within 10 units of each x point. Start by creating some random non-sorted data with N = 10,000 points:

N = 10000; x = randi(300,N,1) + 20+3*randn(N,1) ; y = 3*sind(x) + randn(size(x)) + 3; plot(x,y,'k.') axis tight box off

Now get the moving median of all points within 15 x units of each value:

yb = scatstat1(x,y,15,@median,'omitnan'); hold on plot(x,yb,'bo')

Note, this function does not interpolate or enforce any equal spacing. Final values correspond to each x point, which need not be equally spaced or sorted. Here's a zoom-in so you can see that the x values are not equally spaced:

axis([167 173 1 6])

## Example 2: Sea ice data

From 1978 to 1988 the NSIDC sea ice time series contains about one measurement about every two days. After 1988, it's a daily time series. What if you want a 2 year moving average of sea ice time series, on this dataset which does not have constant temporal resolution?

Here's the situation:

load seaice_extent.mat % loads the sample data figure plot(t,extent_N) axis tight box off ylabel 'northern hemisphere sea ice extent 1\times10^6 km)'

Before 1988, the data is available every other day, and then there's a gap in the data, and then it's daily data. Plot the time between measurements (in days) as the gradient of the `datenum` values:

figure plot(t,gradient(datenum(t))) ntitle 'days between measurements' box off

You see, there's a gap.

For a two-year moving average (a radius of 1 year), convert the datetime format `t` to `datenum` format, and use `scatstat1` with a radius of 365.25 days:

extent_N_2yr = scatstat1(datenum(t),extent_N,365.25); figure(1) % goes back to the sea ice time series figure hold on plot(t,extent_N_2yr,'k','linewidth',2)

And of course, you probably don't want to trust anything within 1 radius (1 year for the 2 year moving average) of the ends, which is why the strange wiggles appear on at the start and end of the time series.

## Author info:

This function was written by Chad A. Greene of the University of Texas at Austin's Institute for Geophysics (UTIG), June 2016.