Find Events in Timestamped Data

This example shows how to use events to work with timestamped data stored in timetables. An event consists of a timestamp (when something happened), often a description (what happened), and sometimes additional information about the event. The timestamp for an instantaneous event is a single point in time, either a datetime or duration value. The timestamp for an interval event consists of two points in time that define the start and end of the event.

You can use events to find or mark data of interest for plotting and analysis. This example shows how you can define events in data using information that is within your data set. To add events from an external data source to your timestamped data, see Add Events from External Data to Timetable.

You can store events separately from the main data set. A timetable is a convenient way to store events because it can hold the event timestamps along with their descriptions and any other information. You can use the event timestamps as subscripts into the main data set to select data at, before, or after an event, or between two events. You can also store events in the same timetable as the data set they describe by using one or more additional timetable variables. This example shows how to use both representations. The representation that is more useful depends on the data analysis you plan to perform.

Define Data Set Using Variations in Length of Day

By definition a day is 86,400 seconds long, where the second has a precise definition in the International System of Units (SI). However, the length of a day actually varies due to several physical causes. It varies with the seasons by as much as 30 seconds over and 21 seconds under the SI definition because of the eccentricity of Earth's orbit and the tilt of its axis. Averaging these seasonal effects enables the definition of the mean solar day, which does not vary in length over a year.

Also, there is a very long-term slowing in the rotational speed of the Earth due to tidal interaction with the moon; a smaller, opposite, shorter-term component believed to be due to melting of continental ice sheets; very short-term cycles on the order of decades; and unpredictable fluctuations due to geological events and other causes. Because of those effects, the length of a mean solar day might increase or decrease. In recent decades, it has fluctuated up and down, but has mostly been 1–3 milliseconds longer than 86,400 seconds. That difference is known as the excess Length of Day, or excess LOD.

For this example, create a timetable that contains the excess LOD for every day from January 1, 1962, to the present. The International Earth Rotation and Reference Systems Service (IERS) collects and publishes this data. However, this data needs preprocessing before storing in a MATLAB timetable because the dates are modified Julian dates. To read the IERS data into a table, use the readtable function. Rename the two variables of interest to MJD and ExcessLOD.

file = "https://datacenter.iers.org/data/latestVersion/223_EOP_C04_14.62-NOW.IAU1980223.txt";
IERSdata = readtable(file,"NumHeaderLines",14);
IERSdata.Properties.VariableNames([4 8]) = ["MJD","ExcessLOD"];

To store the excess LOD values in a timetable, convert the modified Julian dates to datetime values. Use the datetime function with the "ConvertFrom","mjd" name-value argument. Then convert IERSdata from a table to a timetable using the table2timetable function.

IERSdata.Date = datetime(IERSdata.MJD,"ConvertFrom","mjd");
IERSdata.ExcessLOD = seconds(IERSdata.ExcessLOD);
IERSdata = table2timetable(IERSdata(:,["Date","ExcessLOD"]))

IERSdata=22017×1 timetable
       Date         ExcessLOD  
    ___________    ____________

    01-Jan-1962    0.001723 sec
    02-Jan-1962    0.001669 sec
    03-Jan-1962    0.001582 sec
    04-Jan-1962    0.001496 sec
    05-Jan-1962    0.001416 sec
    06-Jan-1962    0.001382 sec
    07-Jan-1962    0.001413 sec
    08-Jan-1962    0.001505 sec
    09-Jan-1962    0.001628 sec
    10-Jan-1962    0.001738 sec
    11-Jan-1962    0.001794 sec
    12-Jan-1962    0.001774 sec
    13-Jan-1962    0.001667 sec
    14-Jan-1962     0.00151 sec
    15-Jan-1962    0.001312 sec
    16-Jan-1962    0.001112 sec
      ⋮

Plot the excess LOD as a function of time.

plot(IERSdata.Date,IERSdata.ExcessLOD,"b-");
ylabel("Excess LOD");

Find Events in Data Set

Since the 1960s there have been several periods when the excess LOD decreased over the short term. If you smooth the excess LOD data, you can see this local behavior more easily.

To smooth the excess LOD, use the smoothdata function. Then plot the smoothed data over the excess LOD.

IERSdata.SmoothedELOD = smoothdata(seconds(IERSdata.ExcessLOD),"loess","SmoothingFactor",.4);
plot(IERSdata.Date,IERSdata.ExcessLOD,"b-");
hold on
plot(IERSdata.Date,IERSdata.SmoothedELOD,"g-","LineWidth",2);
hold off
ylabel("Excess LOD");

The peaks and troughs of the smoothed data show where the short-term trend changed direction. After reaching a peak, the excess LOD decreases. After reaching a trough, the excess LOD increases. The peaks and troughs are notable events in this data set.

To identify the peaks and troughs in the smoothed data, use the islocalmax and islocalmin functions. Then get the date and the value of the excess LOD for each peak and trough. Create a categorical array with two types, peak and trough, that describe these two types of events.

peaks = find(islocalmax(IERSdata.SmoothedELOD));
troughs = find(islocalmin(IERSdata.SmoothedELOD));
Date = IERSdata.Date([peaks;troughs]);
Value = IERSdata.SmoothedELOD([peaks;troughs]);
Type = categorical([zeros(size(peaks)); ones(size(troughs))],[0 1],["peak","trough"]);

Save Events in Timetable

Create a timetable that contains the peak and trough events. Each row of the timetable contains the event date, the type of event, and the value of the excess LOD at that event.

extremaEvents = sortrows(timetable(Date,Type,Value))

extremaEvents=9×2 timetable
       Date         Type       Value   
    ___________    ______    __________

    16-Jul-1972    peak       0.0030145
    04-Aug-1975    trough     0.0027592
    25-Jun-1977    peak       0.0028619
    11-Jan-1987    trough     0.0012623
    31-Oct-1993    peak       0.0023055
    11-Nov-2003    trough    0.00031077
    13-Feb-2008    peak      0.00087789
    18-Jul-2010    trough    0.00074339
    24-Dec-2015    peak       0.0012426

Plot Events Against Data

Mark the peaks and troughs on the plot, using a triangle pointed upward for peaks and a triangle pointed downward for troughs.

hold on
isPeak = (extremaEvents.Type == "peak");
hpeaks = plot(extremaEvents.Date(isPeak),extremaEvents.Value(isPeak),"y^","MarkerFaceColor","y");
isTrough = (extremaEvents.Type == "trough");
htroughs = plot(extremaEvents.Date(isTrough),extremaEvents.Value(isTrough),"yv","MarkerFaceColor","y");
hold off

Find Data Between Two Events

Select the subset of data between the second peak and the second trough.

First specify the time range of these events by using the timerange function and the dates on which the peak and trough occurred.

between2peak2trough = timerange(extremaEvents.Date(3),extremaEvents.Date(4))

between2peak2trough = 
	timetable timerange subscript:

		Select timetable rows with times in the half-open interval:
		[25-Jun-1977 00:00:00, 11-Jan-1987 00:00:00)

	See Select Times in Timetable.

Index into IERSdata using the time range. Display the start and the end of the subset using the head and tail functions.

segment4 = IERSdata(between2peak2trough,:);
head(segment4)

       Date         ExcessLOD      SmoothedELOD
    ___________    ____________    ____________

    25-Jun-1977    0.002277 sec     0.0028619  
    26-Jun-1977    0.002236 sec     0.0028619  
    27-Jun-1977    0.002147 sec     0.0028619  
    28-Jun-1977    0.002032 sec     0.0028619  
    29-Jun-1977    0.001936 sec     0.0028619  
    30-Jun-1977    0.001871 sec     0.0028619  
    01-Jul-1977    0.001875 sec     0.0028619  
    02-Jul-1977     0.00193 sec     0.0028619  

tail(segment4)

       Date          ExcessLOD      SmoothedELOD
    ___________    _____________    ____________

    03-Jan-1987    0.0015672 sec     0.0012623  
    04-Jan-1987     0.001762 sec     0.0012623  
    05-Jan-1987    0.0018557 sec     0.0012623  
    06-Jan-1987    0.0018304 sec     0.0012623  
    07-Jan-1987    0.0016989 sec     0.0012623  
    08-Jan-1987    0.0014924 sec     0.0012623  
    09-Jan-1987    0.0012496 sec     0.0012623  
    10-Jan-1987    0.0010093 sec     0.0012623  

Create Interval Events from Instantaneous Events

From peak to trough, the excess LOD is decreasing, meaning that the Earth's rotation speeds up during that interval. Store these decreasing periods in a table as interval events. The last peak has no matching trough, so use the last time in the data set as the end time for the last interval.

StartDate = IERSdata.Date(peaks);
EndDate = [IERSdata.Date(troughs); IERSdata.Date(end)];
decreasingEvents = table(StartDate,EndDate)

decreasingEvents=5×2 table
     StartDate       EndDate  
    ___________    ___________

    16-Jul-1972    04-Aug-1975
    25-Jun-1977    11-Jan-1987
    31-Oct-1993    11-Nov-2003
    13-Feb-2008    18-Jul-2010
    24-Dec-2015    12-Apr-2022

Compute the average decrease in excess LOD during each interval (in units of seconds of daily excess LOD per year). Add that information to the interval events in the table.

dTime = decreasingEvents.EndDate - decreasingEvents.StartDate;
dExcess = IERSdata.SmoothedELOD(decreasingEvents.EndDate) - IERSdata.SmoothedELOD(decreasingEvents.StartDate);
decreasingEvents.AnnualAvgDecrease = seconds(dExcess ./ years(dTime))

decreasingEvents=5×3 table
     StartDate       EndDate      AnnualAvgDecrease
    ___________    ___________    _________________

    16-Jul-1972    04-Aug-1975     -8.3728e-05 sec 
    25-Jun-1977    11-Jan-1987     -0.00016756 sec 
    31-Oct-1993    11-Nov-2003      -0.0001989 sec 
    13-Feb-2008    18-Jul-2010     -5.5449e-05 sec 
    24-Dec-2015    12-Apr-2022      -0.0002662 sec 

These results show that the mean solar day, averaged over an entire year, has been decreasing over the last few years by about 0.3 milliseconds per year. The mean solar day is currently near or even slightly less than 86,400 seconds. However, many experts believe that this trend will not continue.

Convert Instantaneous Events to State Variable Using retime Function

The extremaEvents timetable records times at which the smoothed excess LOD reached a peak and began a decrease or reached a trough and began an increase. Another way to represent those changes is as a state variable that indicates whether the data are increasing or decreasing at any given time in the data set. To convert the discrete peak and trough events into a new increasing and decreasing state variable, first assign the peak and trough events into the excess LOD data at their time of occurrence.

The new timetable variable, State, is a state variable

with three states. The states are peak at any peak, trough at any trough, and <undefined> at any other time between a peak and a trough.

IERSdata.State(extremaEvents.Date) = extremaEvents.Type;
IERSdata(timerange("2015-12-23","2016-01-01"),:)

ans=9×3 timetable
       Date          ExcessLOD      SmoothedELOD       State   
    ___________    _____________    ____________    ___________

    23-Dec-2015    0.0016556 sec     0.0012426      <undefined>
    24-Dec-2015    0.0015032 sec     0.0012426      peak       
    25-Dec-2015    0.0014242 sec     0.0012426      <undefined>
    26-Dec-2015    0.0014033 sec     0.0012426      <undefined>
    27-Dec-2015    0.0014386 sec     0.0012426      <undefined>
    28-Dec-2015    0.0015334 sec     0.0012426      <undefined>
    29-Dec-2015    0.0016433 sec     0.0012426      <undefined>
    30-Dec-2015    0.0017386 sec     0.0012425      <undefined>
    31-Dec-2015    0.0018316 sec     0.0012425      <undefined>

Now rename the values in State to represent the state at any time, not just at peaks and troughs. Use the renamecats function to rename the states as decreasingLOD and increasingLOD. These two states represent decreasing LOD and increasing LOD.

To carry these state variable values forward to each time in the data, use the retime function, specifying the previous method. After every trough, retime uses the previous method to carry increasingLOD forward as the state, until retime encounters a peak. Then it carries decreasingLOD forward as the state.

IERSdata.State = renamecats(IERSdata.State,["decreasingLOD","increasingLOD"]);
IERSdata = retime(IERSdata,IERSdata.Date,"previous");
IERSdata(timerange("2015-12-23","2016-01-01"),:)

ans=9×3 timetable
       Date          ExcessLOD      SmoothedELOD        State    
    ___________    _____________    ____________    _____________

    23-Dec-2015    0.0016556 sec     0.0012426      increasingLOD
    24-Dec-2015    0.0015032 sec     0.0012426      decreasingLOD
    25-Dec-2015    0.0014242 sec     0.0012426      decreasingLOD
    26-Dec-2015    0.0014033 sec     0.0012426      decreasingLOD
    27-Dec-2015    0.0014386 sec     0.0012426      decreasingLOD
    28-Dec-2015    0.0015334 sec     0.0012426      decreasingLOD
    29-Dec-2015    0.0016433 sec     0.0012426      decreasingLOD
    30-Dec-2015    0.0017386 sec     0.0012425      decreasingLOD
    31-Dec-2015    0.0018316 sec     0.0012425      decreasingLOD

retime is a valuable tool to convert from instantaneous events to a state variable. In this example, the peaks and troughs occurred at certain row times of the data timetable. But even for externally defined events whose times are not in your data set, retime enables you to evaluate the state at only the times that do appear in the data. In this example, passing IERSdata.Date as an argument into retime shows this capability.

Finally, highlight the segments where excess LOD is decreasing in red. Even though decreasingEvents contains the start and end times of those segments, you can plot those segments more easily using the state variable.

delete([hpeaks; htroughs]);
hold on
decreasing = (IERSdata.State == "decreasingLOD");
plot(IERSdata.Date(decreasing),IERSdata.SmoothedELOD(decreasing),'r.');
hold off

Alternatively, highlight the background in those regions. In this case, the interval events are more convenient.

hold on
startEnd = [decreasingEvents.StartDate decreasingEvents.EndDate];
h = fill(startEnd(:,[1 2 2 1]),[-.002 -.002 .005 .005],"red","FaceAlpha",.2,"LineStyle","none");
hold off

Find More Complex Events in Data

The excess LOD has both increased and decreased since the 1960s. Indeed, in many years there were short periods when the raw excess LOD was significantly negative. These are only very short-term fluctuations, but during those periods the Earth was rotating one millisecond or more faster than 86,400 SI seconds.

plot(IERSdata.Date,IERSdata.ExcessLOD,"b-");
ylabel("Excess LOD");
hold on
line(IERSdata.Date([1 end]),[0 0],"Color","k","LineStyle",":")
hold off
ylabel("Excess LOD");

Identify the years in which the excess LOD was negative on any day. Then use retime to find the minimum excess LOD in each of those years and create a timetable of events. These are interval events in one sense but are stored as instantaneous events marked only by their year.

negLOD = IERSdata(IERSdata.ExcessLOD < 0,"ExcessLOD");
negYears = unique(dateshift(negLOD.Date,"start","year"));
negYears.Format = "uuuu";
negLODEvents = retime(negLOD,negYears,"min");
negLODEvents.Properties.VariableNames = "MinExcessLOD"

negLODEvents=25×1 timetable
    Date     MinExcessLOD 
    ____    ______________

    1984     -9.38e-05 sec
    1986     -1.33e-05 sec
    1987    -0.0001492 sec
    1988     -7.06e-05 sec
    1999    -0.0001063 sec
    2000     -0.000311 sec
    2001    -0.0007064 sec
    2002    -0.0007436 sec
    2003    -0.0009769 sec
    2004    -0.0010672 sec
    2005    -0.0010809 sec
    2006    -0.0003865 sec
    2007    -0.0006192 sec
    2008    -0.0003945 sec
    2009    -0.0004417 sec
    2010     -0.000784 sec
      ⋮

Mark the time axis red for each year that had periods when the excess LOD was negative. In this data set, such years happen more frequently after the year 2000.

hold on
plot([negLODEvents.Date negLODEvents.Date+calyears(1)],[-.0016 -.0016],"r-","lineWidth",6);
ylim(seconds([-.0016 .0045]));
hold off

Whether you use instantaneous events, state variables, or switch between them depends on which representation is more convenient and useful for the data analysis that you plan to carry out. Both representations are useful ways to add information about events to your timestamped data in a timetable.