MODIS
M.Moody
2025-11-06
Introduction MODIS data refers to the information collected by the Moderate Resolution Imaging Spectroradiometer (MODIS) instrument aboard NASA’s two different satellites:Terra and Aqua. Terra’s orbit around the Earth is timed so that it passes from north to south across the equator in the morning, while Aqua passes south to north over the equator in the afternoon. Terra MODIS and Aqua MODIS are viewing the entire Earth’s surface every 1 to 2 days, acquiring data in 36 spectral bands!
MODIS Vegetation Data MODIS vegetation indices, are produced on 16-day intervals and at multiple spatial resolutions. They provide consistent spatial and temporal comparisons of vegetation canopy greenness, a composite property of leaf area, chlorophyll and canopy structure. Two vegetation indices are derived from atmospherically-corrected reflectance in the red, near-infrared, and blue wavebands; the normalized difference vegetation index (NDVI), and the enhanced vegetation index (EVI), which minimizes canopy-soil variations and improves sensitivity over dense vegetation conditions.
We are going to create a MODIS NDVI time series for your research locations. you can obtain MODIS NDVI data from sources like Google Earth Engine or NASA Earthdata, process it by defining your area of interest, filtering and grouping the data, and performing temporal smoothing to reduce noise. You can then analyze the smoothed time series to study vegetation phenology, productivity, or changes over time by identifying trends or anomalies.
We could use GEE or some other repository to get the MODIS data, but today we will use the MODISTools package in R with subsets data from the Oak Ridge National Laboratory Distributed Active Archive Center. This allows us to obtain MODIS time series for specific locations or small regions of interest directly within R.
Before beginning set your working directory. And put all files you will need in that directory.
Exercise on Modis data
Start by loading packages.
We can then explore MODIS products.
products <- mt_products() # All products available
head(products) #Look at product 43## product
## 1 Daymet
## 2 ECO4ESIPTJPL
## 3 ECO4WUE
## 4 GEDI03
## 5 GEDI04_B
## 6 MCD12Q1
## description
## 1 Daily Surface Weather Data (Daymet) on a 1-km Grid for North America, Version 4 R1
## 2 ECOSTRESS Evaporative Stress Index PT-JPL (ESI) Daily L4 Global 70 m
## 3 ECOSTRESS Water Use Efficiency (WUE) Daily L4 Global 70 m
## 4 GEDI Gridded Land Surface Metrics (LSM) L3 1km EASE-Grid, Version 2
## 5 GEDI Gridded Aboveground Biomass Density (AGBD) L4B 1km EASE-Grid, Version 2.1
## 6 MODIS/Terra+Aqua Land Cover Type (LC) Yearly L3 Global 500 m SIN Grid
## frequency resolution_meters
## 1 1 day 1000
## 2 Varies 70
## 3 Varies 70
## 4 One time 1000
## 5 One time 1000
## 6 1 year 500Let’s get a listing of the bands of the vegetation indices product. Product 43 in the list is VNP13A1. VNP13A1 is a dataset that provides 16-day global vegetation index information at 500-meter resolution, derived from the VIIRS sensor on the Suomi-NPP satellite. This dataset includes (NDVI) and (EVI).
bands <- mt_bands(product = "VNP13A1")
head(bands)## band
## 1 500_m_16_days_blue_reflectance
## 2 500_m_16_days_composite_day_of_the_year
## 3 500_m_16_days_EVI
## 4 500_m_16_days_EVI2
## 5 500_m_16_days_green_reflectance
## 6 500_m_16_days_NDVI
## description units valid_range
## 1 Blue band (M3 478-498 nm) reflectance reflectance 0 to 10000
## 2 Day of the year Julian day of the year 1 to 366
## 3 16 day EVI average EVI ratio - No units -10000 to 10000
## 4 16 day EVI2 average EVI2 ratio - No units -10000 to 10000
## 5 Green band (M4 545-565 nm) reflectance reflectance 0 to 10000
## 6 16 day NDVI average NDVI ratio - No units -10000 to 10000
## fill_value scale_factor add_offset
## 1 -1000 0.0001 0
## 2 -1 <NA> <NA>
## 3 -15000 0.0001 0
## 4 -15000 0.0001 0
## 5 -1000 0.0001 0
## 6 -15000 0.0001 0Define the location and examine dates for the remote sensing product.
dates <- mt_dates(product = "VNP13A1", lat = 39.0542734, lon = -79.6700617)
head(dates)## modis_date calendar_date
## 1 A2012017 2012-01-17
## 2 A2012025 2012-01-25
## 3 A2012033 2012-02-02
## 4 A2012041 2012-02-10
## 5 A2012049 2012-02-18
## 6 A2012057 2012-02-26tail(dates)## modis_date calendar_date
## 564 A2024105 2024-04-14
## 565 A2024113 2024-04-22
## 566 A2024121 2024-04-30
## 567 A2024129 2024-05-08
## 568 A2024137 2024-05-16
## 569 A2024145 2024-05-24Download NDVI data for White Mountains
How long does it take to download 200 km2 of VIIRS NDVI data at 500 m resolution for 2 bands?
end_time - start_time## Time difference of 28.37962 secsFinish up by summarizing and plotting the data.
WM_ndvi %>%
filter(band == "500_m_16_days_NDVI") %>%
group_by(calendar_date) %>%
summarise(doy = yday(as_date(calendar_date)),
ndvi_median = median(value * as.numeric(scale))) %>%
distinct(doy, .keep_all = TRUE) -> WM_med_ndvi## Warning: Returning more (or less) than 1 row per `summarise()` group was deprecated in
## dplyr 1.1.0.
## ℹ Please use `reframe()` instead.
## ℹ When switching from `summarise()` to `reframe()`, remember that `reframe()`
## always returns an ungrouped data frame and adjust accordingly.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## `summarise()` has grouped output by 'calendar_date'. You can override using the
## `.groups` argument.head(WM_med_ndvi)## # A tibble: 6 × 3
## # Groups: calendar_date [6]
## calendar_date doy ndvi_median
## <chr> <dbl> <dbl>
## 1 2023-01-01 1 0.414
## 2 2023-01-09 9 0.350
## 3 2023-01-17 17 0.490
## 4 2023-01-25 25 0.444
## 5 2023-02-02 33 0.489
## 6 2023-02-10 41 0.488Use GGPlot to make the Plot
