U.S. Department of the Interior

U.S. Geological Survey

Open-File Report 2016–1120

A Satellite Model of Southwestern Willow Flycatcher 

(Empidonax traillii extimus) Breeding Habitat and a 

Simulation of Potential Effects of Tamarisk Leaf Beetles 

(Diorhabda spp.), Southwestern United States

A Satellite Model of Southwestern Willow Flycatcher

(Empidonax traillii extimus) Breeding Habitat and a 

Simulation of Potential Effects of Tamarisk Leaf Beetles

(Diorhabda spp.), Southwestern United States

By James R. Hatten

Open-File Report 2016–1120 

U.S. Department of the Interior

U.S. Geological Survey
U.S. Department of the Interior

SALLY JEWELL, Secretary

U.S. Geological Survey

Suzette M. Kimball, Director

U.S. Geological Survey, Reston, Virginia: 2016

For more information on the USGS—the Federal source for science about the Earth,

its natural and living resources, natural hazards, and the environment—visit

http://www.usgs.gov/ or call 1–888–ASK–USGS (1–888–275–8747).

For an overview of USGS information products, including maps, imagery, and publications,

visit http://www.http://www.store.usgs.gov.

Any use of trade, firm, or product names is for descriptive purposes only and does not imply

endorsement by the U.S. Government.

Although this information product, for the most part, is in the public domain, it also may

contain copyrighted materials as noted in the text. Permission to reproduce copyrighted items

must be secured from the copyright owner.

Suggested citation:

Hatten, J.R., 2016, A satellite model of Southwestern Willow Flycatcher (Empidonax traillii extimus) breeding 

habitat and a simulation of potential effects of tamarisk leaf beetles (Diorhabda spp.), Southwestern United States: 

U.S. Geological Survey Open-File Report 2016–1120, 88 p., http://dx.doi.org/10.3133/ofr20161120. 

ISSN 2331-1258 (online)] 
iii

Contents

Executive Summary.................................................................................................................................................... 1

Introduction................................................................................................................................................................. 2

Methods...................................................................................................................................................................... 6

Satellite-Based Habitat Predictions......................................................................................................................... 6

Geospatial Database........................................................................................................................................... 6

Application of the Satellite Model ...................................................................................................................... 12

Accuracy Assessment ....................................................................................................................................... 13

Summarization and Display of Model Results ................................................................................................... 14

Habitat Change Detection ................................................................................................................................. 14

Tamarisk Leaf Beetle Impact Assessment ............................................................................................................ 14

Habitat Time Series........................................................................................................................................... 15

Virgin River Beetle-Impact Analysis................................................................................................................... 16

Lower Colorado River Beetle-Impact Simulation ............................................................................................... 18

Upper Gila River Beetle-Impact Simulation ....................................................................................................... 20

Results...................................................................................................................................................................... 22

Satellite-Based Habitat Predictions....................................................................................................................... 22

Rangewide ........................................................................................................................................................ 22

U.S. Fish and Wildlife Service Management Units ............................................................................................ 23

7.5-Minute Quadrangles.................................................................................................................................... 27

U.S. Fish and Wildlife Service Critical-Habitat Reaches.................................................................................... 32

Model Accuracy................................................................................................................................................. 35

Habitat Change Detection ................................................................................................................................. 47

Tamarisk Leaf Beetle Impact Assessment ............................................................................................................ 51

Lower Virgin River............................................................................................................................................. 51

Upper Gila River................................................................................................................................................ 61

Discussion ................................................................................................................................................................ 65

Key Findings and a Path Forward......................................................................................................................... 65

Habitat Change Detection ..................................................................................................................................... 66

Rangewide Ranking of Flycatcher Breeding Habitat............................................................................................. 68

Accuracy Assessment........................................................................................................................................... 68

Beetle-Impact Analysis and Simulations ............................................................................................................... 69

Next Steps ............................................................................................................................................................ 70

Occupancy Modeling......................................................................................................................................... 70

Automated Cloud-Based Mapping..................................................................................................................... 70

Acknowledgments..................................................................................................................................................... 71

References Cited...................................................................................................................................................... 71

Appendix A. Metadata for Landsat Scenes Used in Regionwide Mapping and Habitat Time Series ....................... 75

Appendix B. Normalized Difference Vegetation Index (NDVI) Conversions for Landsat 8....................................... 83

Figures

Figure 1. Maps showing location of the project area and the 35 U.S. Fish and Wildlife Service Southwestern 

Willow Flycatcher management unit boundaries, designated critical-habitat reaches, and waterbodies

considered by the satellite model for regionwide modeling, Southwestern United States.......................................... 3
iv

Figure 2. Map showing yearly distribution of tamarisk beetles (Diorhabda spp.) throughout 

Western United States, 2007–14. (Map produced by the Tamarisk Coalition.).......................................................... 5

Figure 3. Map showing paths and rows of 57 Landsat scenes overlapping breeding range of

Southwestern Willow Flycatcher, Southwestern United States.................................................................................. 7

Figure 4. Map showing all fourth-order or larger streams and waterbodies, plus lakes greater than 1 square 

kilometer, and two elevation zones (1) below 1,524 m, and (2) above 1,524 m. ....................................................... 9

Figure 5. Map showing urban and agriculture areas identified and masked (that is, excluded) from the 

satellite model.......................................................................................................................................................... 10

Figure 6. Maps showing effects of masking on output from satellite model at confluence of the Colorado 

and Gila Rivers, Arizona and California ................................................................................................................... 11

Figure 7. Map showing locations of 758 Southwestern Willow Flycatcher (SWFL) territories that were

used for model verification overlaid on 35 U.S. Fish and Wildlife Service management unit boundaries, 

Southwestern United States, summer 2014............................................................................................................. 13

Figure 8. Map showing boundary (polygon with seven segments) of lower Virgin River where a 

habitat time series was calculated and a change detection was done, Nevada and Arizona, 1986–2015............... 17

Figure 9. Map showing Lower Colorado River boundary, tamarisk patches and structure (I–VI), U.S. Fish and 

Wildlife Service (FWS) Southwestern Willow Flycatcher (SWFL) management unit boundaries overlaid on a 

normalized difference vegetation index (NDVI) image obtained July 2015, California, Nevada, and Arizona.......... 19

Figure 10. Map showing project boundary, nine river reaches used in tamarisk leaf beetle-impact simulation, 

and nine vegetation classes, upper Gila River, Arizona, 1986–2015....................................................................... 21

Figure 11. Graph showing area of predicted flycatcher breeding habitat below 1,524-meter elevation and

for all elevations by State. Number above each bar is area of predicted habitat in square kilometers.................... 22

Figure 12. Graphs showing area of predicted flycatcher breeding habitat in U.S. Fish and Wildlife Service 

management units below and above 1,524-meter elevation, Southwestern United States, 2013–15...................... 24

Figure 13. Graph showing changes in predicted flycatcher breeding habitat in U.S. Fish and Wildlife Service 

management units below and above 1,524-meter elevation and for all elevations, 2013–15 .................................. 26

Figure 14. Map showing area of predicted flycatcher breeding habitat (averaged across 2013–15) at 

all elevations in 6,521 U.S. Geological Survey 7.5-minute quadrangles, Southwestern United States.................... 28

Figure 15. Graph showing changes in predicted flycatcher habitat for the 50 top-ranked 

U.S. Geological Survey 7.5-minute quadrangles below 1,524-meter elevation, 2014–15........................................ 29

Figure 16. Maps showing predicted flycatcher breeding habitat at (A) Paraje Well and San Marcial, 

New Mexico, June 2013, and (B) San Carlos Reservoir, south-central Arizona ...................................................... 30

Figure 17. Graph showing density of flycatcher territories (n = 747) in 2014 in five probability 

classes output by the satellite model. ...................................................................................................................... 36

Figure 18. Flycatcher habitat suitability maps of the upper Owens River, California............................................... 37

Figure 19. Flycatcher habitat suitability maps of Alamo Lake, Arizona ................................................................... 39

Figure 20. Flycatcher habitat suitability maps of Rio Grande (Elephant Butte Reservoir), New Mexico.................. 41

Figure 21. Fflycatcher habitat suitability maps of Havasu National Wildlife Refuge (Topock Marsh), Arizona 

and California........................................................................................................................................................... 43

Figure 22. Flycatcher habitat suitability maps of the Santa Margarita and San Luis Rey Rivers

(Camp Pendleton), California................................................................................................................................... 45

Figure 23. Map showing relative changes in predicted flycatcher habitat, as determined from a

satellite model at a 40-percent probability threshold, Southwestern United States, 2013–15.................................. 48

Figure 24. Map showing changes in predicted flycatcher habitat at Lake McMillan, Texas, 

June 2014–June 2015 ............................................................................................................................................. 49

Figure 25. Map showing changes in predicted flycatcher habitat in Fort Craig and San Marcial quadrangles, 

New Mexico, 2014–15 ............................................................................................................................................. 50
v 

Figure 26. A habitat time series of the lower Virgin River (see fig. 8), Nevada and Arizona, produced by running 

the satellite model repeatedly from 1986 to 2015 (2012 omitted) ............................................................................ 53

Figure 27. Map showing changes in predicted flycatcher habitat along the lower Virgin River, Nevada

and Arizona, 2010–2015, as determined from a satellite model at a 40-percent probability threshold. ................... 54

Figure 28. Graph showing changes in area of predicted flycatcher habitat in the lower Virgin River, 

Nevada and Arizona (fig. 27), 2010–2015. .............................................................................................................. 55

Figure 29. Graph showing changes in NDVI values in four zones from 2010 (June 22) to 2015 (June 25—) 

in the lower Virgin River, Nevada and Arizona......................................................................................................... 55

Figure 30. Graph showing changes in predicted flycatcher habitat inside seven zones (fig. 27) along the lower 

Virgin River, Nevada and Arizona, 2010–2015........................................................................................................ 56

Figure 31. Graph showing simulation results showing how tamarisk leaf beetles may affect flycatcher habitat at 

northern part (Topock Marsh) of Havasu National Wildlife Refuge along lower Colorado River, Arizona and 

California ................................................................................................................................................................. 58

Figure 32. Map showing simulation results showing how tamarisk leaf beetles may affect flycatcher habitat at 

Cibola National Wildlife Refuge along lower Colorado River, Arizona and California .............................................. 59

Figure 33. Graph showing three states of flycatcher habitat related to tamarisk leaf beetle-impact 

simulation along lower Colorado River, Arizona and California ............................................................................... 60

Figure 34. Graph showing changes in predicted flycatcher habitat along the lower Colorado River, 

by U.S. Fish and Wildlife Service management units(see fig. 9), after simulating a tamarisk leaf beetle

infestation, Arizona and California ........................................................................................................................... 60

Figure 35. A habitat time series of predicted flycatcher habitat for upper Gila River, Arizona, from 1986

to 2015, plus the area of predicted flycatcher habitat after simulating a tamarisk leaf beetle infestation in

2015 (2015b) .......................................................................................................................................................... 62

Figure 36. Map showing changes in predicted flycatcher habitat along the upper Gila River after conducting a 

tamarisk leaf beetle-impact simulation..................................................................................................................... 63

Figure 37. Graph showing three states of flycatcher habitat produced by a tamarisk leaf beetle-impact 

simulation along the upper Gila River...................................................................................................................... 64

Figure 38. Graph showing simulated changes (percent) in predicted flycatcher habitat within 9 zones 

along the upper Gila River....................................................................................................................................... 64

Figure 39. Drought maps portray a drought that intensifies to the west and decreases in the eastern 

portion of the study area from 2013 to 2015 (U.S. Department of Agriculture, 2015). ............................................. 67

Tables

Table 1. Changes in predicted flycatcher habitat in U.S. Fish and Wildlife Service management units, 

Southwestern United States, 2013–15 .................................................................................................................... 25

Table 2. Top-ranked 50 7.5-minute U.S. Geological Survey quadrangles by average area (hectares) 

of predicted flycatcher habitat, Southwestern United States, 2013–15.................................................................... 31

Table 3. Predicted flycatcher habitat in 87 U.S. Fish and Wildlife Service critical-habitat reaches (Unit name), 

located within 6 recovery units and 35 management units, Southwestern United States, 2013–15 ........................ 33

Table 4. Area of predicted flycatcher habitat before and after tamarisk leaf beetle-impact analysis 

was done for the lower Virgin River, Nevada and Arizona, 2010 and 2015 ............................................................. 52

Table 5. Changes in predicted Southwestern Willow Flycatcher habitat along lower Colorado River, 

California, Nevada, and Arizona, after a beetle-impact simulation was done in 2015.............................................. 57

Table 6. Area of predicted flycatcher habitat before and after a tamarisk leaf beetle-impact simulation 

was done for the upper Gila River, Arizona. ............................................................................................................ 65
vi 

Conversion Factors

Inch/Pound to International System of Units

Multiply By To obtain

Area

acre 0.4047 hectare (ha)

acre 0.4047 square hectometer (hm2) 

acre 0.004047 square kilometer (km2)

International System of Units to Inch/Pound

Multiply By To obtain

Length

meter (m) 3.281 foot (ft) 

kilometer (km) 0.6214 mile (mi)

Area

square meter (m2) 0.0002471 acre 

square meter (m2) 10.76 square foot (ft2) 

hectare (ha) 2.471 acre

hectare (ha) 0.003861 square mile (mi2) 

square hectometer (hm2) 2.471 acre

square kilometer (km2) 247.1 acre

square kilometer (km2) 0.3861 square mile (mi2)

Datums 

Vertical coordinate information is referenced to the North American Vertical Datum of 1988 (NAVD 88).

Horizontal coordinate information is referenced to the North American Datum of 1983 (NAD 83).

Elevation, as used in this report, refers to distance above the vertical datum.

Abbreviations

CI confidence interval

FWS U.S. Fish and Wildlife Service

GAP Gap Analysis Program

GIS geographic information system

NDVI normalized difference vegetation index

Reclamation Bureau of Reclamation

TOA top-of-atmosphere
1 

A Satellite Model of Southwestern Willow Flycatcher 

(Empidonax traillii extimus) Breeding Habitat and a 

Simulation of Potential Effects of Tamarisk Leaf Beetles

(Diorhabda spp.), Southwestern United States

By James R. Hatten

Executive Summary

The study described in this report represents the first time that a satellite model has been used to 

identify potential Southwestern Willow Flycatcher (Empidonax traillii extimus) (hereinafter referred to 

as “flycatcher”) breeding habitat rangewide for 2013–15. Fifty-seven Landsat scenes were required to 

map the entire range of the flycatcher, encompassing parts of six States and more than 1 billion 30-

meter pixels. Predicted flycatcher habitat was summarized in a hierarchical fashion from largest to 

smallest: regionwide, State, U.S. Fish and Wildlife Service (FWS) management unit, 7.5-minute 

quadrangle, and critical-habitat reach. The term “predicted habitat” is used throughout this report to 

distinguish areas the satellite model predicts as suitable flycatcher habitat from what may actually exist 

on the ground. A rangewide accuracy assessment was done with 758 territories collected in 2014, and 

change detection was done with yearly habitat maps to identify how and where habitat changed over 

time. Additionally, effects of tamarisk leaf beetles (Diorhabda spp.) on flycatcher habitat were 

summarized for the lower Virgin River from 2010 to 2015, and simulations of how tamarisk leaf beetles 

may affect flycatcher habitat in the lower Colorado and upper Gila Rivers were done for 2015. Model 

results indicated that the largest areas of predicted flycatcher habitat at elevations below 1,524 meters 

were in New Mexico and Arizona, areas followed in descending order by California, Texas, Nevada, 

Utah, and Colorado. By FWS management unit, the largest area of flycatcher habitat during all 3 years

were the Middle Rio Grande (New Mexico), followed by the Upper Gila (Arizona and New Mexico)

and Middle Gila/San Pedro (Arizona) management units. The area of predicted flycatcher habitat varied 

considerably in 7.5-minute quadrangles, ranging from 0 to1,398 hectares (ha). Averaged across 3 years, 

the top three producing quadrangles were Paraje Well (New Mexico), San Marcial (New Mexico), and 

San Carlos Reservoir (Arizona). The top three FWS critical-habitat reaches in 2015 were Rio Grandemiddle (9,544 ha), San Pedro River (1,779 ha), and Gila River-mid San Carlos (1,356 ha); this ranking 

did not change in 2013 or 2014. Change detection among years showed a large shift in predicted 

flycatcher habitat influenced by drought patterns, with California habitat decreasing and New Mexico 

habitat increasing. An accuracy assessment indicated that 88 percent of territories were correctly 

classified at a 40 percent probability threshold, with an exponential relationship between territory 

densities and five probability classes. A spatially explicit analysis indicated that beetles decreased 

predicted flycatcher habitat 94.2 percent from 2010 to 2015 along the lower Virgin River, with only 5.8 

percent persisting. In contrast, beetle simulations indicated that 64.1 percent of habitat will persist along 

the lower Colorado River and 45 percent will persist along the upper Gila River. This project shows that 

the satellite model adequately predicts flycatcher habitat rangewide, but it lacks the ability to predict 
2 

which patches will be occupied in a given year. The next logical step is the development of an 

occupancy model that ties the habitat predictions of the satellite model to patch occupancy so managers 

can better allocate their resources for survey and restoration activities. Finally, the methods presented in 

this report seem well suited for automated mapping applications and cloud-based resources.

Introduction

The Southwestern Willow Flycatcher (Empidonax traillii extimus, hereinafter referred to as 

“flycatcher”) is a federally endangered species that occurs patchily along rivers and streams in the 

Southwestern United States during May–September (U.S. Fish and Wildlife Service, 1995, 2002). The 

flycatcher breeding habitat is characterized by a mosaic of relatively dense tree and shrub growth, 

typically in association with surface water or saturated soil, interspersed with more open areas, open 

water, or shorter, sparser vegetation along rivers, streams, or other wetlands. Plant species composition, 

vegetation height and density, and patch size vary greatly, but most occupied sites typically consist of 

dense vegetation in the interior of the patch and within 3–4 m of the ground (Sogge and Marshall, 2000; 

U.S. Fish and Wildlife Service, 2002; Sogge and others, 2010). Because habitat loss and degradation are 

the main factors contributing to the decline of the species, the flycatcher recovery plan emphasizes the 

increase and improvement of breeding habitat through restoration of native breeding habitat and the 

management of exotic vegetation (U.S. Fish and Wildlife Service, 2002). The estimated rangewide 

breeding population of flycatchers as of the last compilation in 2007 was roughly 1,299 territories 

distributed among 288 known sites (Durst and others, 2008). 

In 2013, the U.S. Fish and Wildlife Service (FWS) designated revised critical habitat for the 

flycatcher under the Endangered Species Act (U.S. Fish and Wildlife Service, 2013). A total of about 

1,975 stream kilometers were designated as critical habitat, located in 6 states, 38 counties, 24 

management units, and 6 recovery units. Identified by stream segment, the lateral extent includes 

riparian areas and streams in the 100-year floodplain, encompassing about 84,569 ha (fig. 1). The FWS 

would like to acquire the ability to locate, map, and monitor potential flycatcher breeding habitat 

rangewide, including areas not designated as critical habitat. Most flycatcher habitat is inherently 

dynamic, with individual riparian patches subject to cycles of creation, growth, and loss because of

drought, flooding, fire, and other disturbances. Former breeding patches can lose suitability quickly, 

whereas new habitat can develop in a few years, especially in reservoir drawdown zones (Paradzick and 

Hatten, 2004; Paxton and others, 2007). Because the distribution and extent of flycatcher habitat 

changes over time, an accurate understanding and assessment of current habitat suitability requires upto-date information on vegetation characteristics and other variables associated with occupied sites. 
3 

Figure 1. Maps showing location of the project area and the 35 U.S. Fish and Wildlife Service Southwestern 

Willow Flycatcher management unit boundaries, designated critical-habitat reaches, and waterbodies considered 

by the satellite model for regionwide modeling, Southwestern United States.
4 

The Arizona Game and Fish Department developed a flycatcher remote sensing/geographic 

information system (GIS) habitat model (hereinafter referred to as the “satellite model”) that uses 

Landsat Thematic Mapper imagery and a 30-m-resolution digital elevation model (Hatten and 

Paradzick, 2003). The satellite model was developed with presence/absence survey data acquired along 

the San Pedro and Gila Rivers, and from Salt River and Tonto Creek inlets to Roosevelt Lake in 

southern Arizona. The satellite model uses a logistic regression equation to divide riparian vegetation 

into a continuous range of probabilities extending from almost 0 to 99 percent, with higher probabilities 

most likely to contain a flycatcher territory. The model has been successfully tested at Alamo Lake, 

Arizona (Hatten and Paradzick, 2003); applied statewide in Arizona with 2001 imagery (Dockens and 

others, 2004); and projected along the Rio Grande River, New Mexico (Hatten and Sogge, 2007). In 

each case, the satellite model performed within expectations by identifying riparian areas with the 

highest densities of flycatcher territories. In a follow-up study, the satellite model explained 79 percent

of the fluctuation in the flycatcher breeding population at Roosevelt Lake from 1996 to 2005 (Paxton 

and others, 2007; Hatten and others, 2010). 

A secondary objective of this project was to identify how and where the recently introduced 

tamarisk leaf beetle (Diorhabda spp.) is affecting flycatcher habitat in the Southwestern United States. 

In 2001, tamarisk leaf beetles (hereinafter referred to as “beetle”) were released by the U.S. Department 

of Agriculture Animal and Plant Health Inspection Service at 10 sites in six States (California, Nevada, 

Utah, Colorado, Wyoming, and Texas) to control invasive tamarisk (Tamarix spp.) (Bean and others, 

2013). Tamarisk, which is an introduced phreatophyte and a common nesting substrate for flycatchers 

and other birds (Sogge and others, 2008; Paxton and others, 2011), spread rapidly following its release 

over a century ago and is now a dominant riparian species throughout the West (Chew, 2013). Feeding 

by the beetle larvae defoliates the tamarisk during the growing season and reduces plant vigor, 

sometimes resulting in plant mortality within 5 years (Hultine and others, 2009, 2015). Since their 

introduction, beetles have spread naturally into most drainages of the upper Colorado River Basin 

(Jamison and others, 2015), southern Utah and Nevada, and the main stem lower Colorado River in 

Arizona (fig. 2). 
5 

Figure 2. Map showing yearly distribution of tamarisk beetles (Diorhabda spp.) throughout Western United States, 

2007–14. (Map produced by the Tamarisk Coalition.)

The northern tamarisk leaf beetle (D. carinulata) has colonized sites in the northern part of the 

Rio Grande and the subtropical tamarisk beetle (D. sublineata) is spreading from western Texas, 

colonizing the southern part of the Rio Grande from Hatch, New Mexico, south to the Mexico border 

(Johnson and Jamison, 2015). The beetle species originating from western Texas are known to have 

additional generations (4–5 generations), leading to an increase in distribution and defoliation as 

population growth expands throughout the Rio Grande corridor. The beetle is expected to be prevalent 

throughout the Lower Colorado, Rio Grande, and Gila River systems within the next 2–3 years ( J. 

Tracy, Texas A&M University, written commun., April 2016), where tamarisk is widely used for 

breeding by flycatchers and other riparian birds. The pace and extent of beetle expansion may cause 

widespread damage to flycatcher breeding habitat over the coming decades (Paxton and others, 2011; 

Bean and others, 2013; Sogge and others, 2013), as it already has caused along the lower Virgin River 

(Bateman and others, 2013; Nagler and others, 2014). The ability to simulate potential effects of beetles 

on flycatcher habitat years in advance would provide managers with a head start in restoring and 

enhancing river reaches before their arrival, thereby minimizing damage to critical flycatcher habitat. 
6 

In the face of so much dynamism, the only effective way to locate, map, and monitor flycatcher 

habitat rangewide is with a combination of remote sensing, field surveys, and GIS technology (U.S. Fish 

and Wildlife Service, 2002; Hatten and others, 2010). The satellite model has broad application for 

identifying and classifying flycatcher breeding habitat rangewide. However, details of flycatcher habitat 

characteristics can vary in different parts of the range, especially between low- and high-elevation sites. 

Additionally, the distinctness of riparian patches and adjacent natural and human-influenced habitat 

differs across the range. Therefore, applying the satellite model rangewide requires significant 

development of a geospatial database, testing and refinement of techniques used to automate the 

application and processing of satellite imagery, and model verification in different parts of flycatcher 

range. This project had three distinct objectives: 

1. To identify, map, and monitor flycatcher breeding habitat rangewide with a published 

satellite model (Hatten and Paradzick, 2003); 

2. To identify, quantify, and monitor beetle damage to flycatcher breeding habitat with the 

satellite model; and

3. To develop a methodology to simulate beetle impacts to flycatcher habitat so managers can 

identify potential restoration and enhancement strategies in advance. 

Methods

Satellite-Based Habitat Predictions

Geospatial Database

The satellite model requires 57 Landsat scenes in order to map “predicted flycatcher breeding 

habitat” rangewide (fig. 3), comprising about 1.1 billion 900-m2 pixels (1,002,839 km2). The word

“predicted” is used throughout this report to distinguish areas the satellite model predicts as suitable 

flycatcher habitat from what areas may actually exist on the ground. Additionally, 30 Landsat scenes 

were used (1986–2015) when producing habitat time series at numerous locations. Producing a 

rangewide map and developing habitat time series required careful attention to image preparation and 

numerical adjustments because of differences in Landsat sensors and acquisition dates (Hatten and 

others, 2010; Hatten, 2014), which was accomplished in several steps. First, differences in predicted 

habitat owing to seasonal effects were minimized by selecting images that were acquired as close to 

June 21 as possible. It was not possible to select the same date each year because clouds sometimes 

obscured the view, or image quality sometimes varied because of satellite malfunctions. In such 

circumstances, the next available Landsat scene was selected, which was 16 days prior to or after the 

target date, because Landsat reimages the Earth every 16 days. Each image was previewed before 

selection with the U.S. Geological Survey GLOVIS tool (U.S. Geological Survey, 2005). The same 

sensor (Landsat 5 Thematic Mapper) was operational from 1986 to 2011, so those scenes were entered 

directly into the satellite model without applying any adjustments to the imagery (see 

appendix A (A1–A3) for Landsat imagery metadata related to rangewide modeling from 2013 to 2015). 
7 

When using Landsat 8 imagery, adjustments were made for a numerical scaling issue (Landsat 5 is 8 bit 

and Landsat 8 is 16 bit) and for differences in normalized difference vegetation index (NDVI) values 

because of slight differences in the red and reflected infrared bands. Specifically, Landsat 8 scenes were 

adjusted (years 2013–15) to top-of-atmosphere (TOA) reflectance to correct the numerical scaling issue 

between the two sensors (Chander and others, 2009). Landsat 8 NDVI values then were adjusted with a 

two-step linear regression so they matched Landsat 5 values; without these two adjustments, the satellite 

model will not work correctly with Landsat 8 imagery (see Chander and others, 2009, for TOA 

corrections and appendix B for details about the two-step regression for NDVI).

Figure 3. Map showing paths and rows of 57 Landsat scenes overlapping breeding range of Southwestern Willow 

Flycatcher, Southwestern United States.
8 

The scale of the modeling required development of new techniques to identify potential 

flycatcher habitat and to create the most specific and accurate maps possible. 

1. Riparian areas where flycatcher territories may occur were identified by buffering streams, 

rivers, and waterbodies according to their size, 

2. All urban and agricultural areas were masked (excluded) because flycatchers do not breed in 

such areas, 

3. Anthropogenic features (for example, circles, rectangles, lines) were identified and masked on 

Landsat scenes that were missed with the distance-to-water mask, and

4. All potentially suitable areas above or below 1,524-m elevation were identified and attributed

because this is an important elevation boundary (Hatten and Paradzick, 2003). 

A variable-width buffer was created for modeling purposes in several steps. First, the 

river/stream reaches and waterbodies nearest the 288 known flycatcher sites were identified by 

overlaying flycatcher sites and the NHDPlus2 stream network (McKay and others, 2012). This 

exploratory analysis indicated that most flycatcher sites were adjacent to stream or river reaches greater

than third order (fourth to ninth order) and lakes greater than 1 km2. Second, given that larger-order 

streams or rivers have larger floodplains, and thus more riparian area that may contain a flycatcher 

territory, a variable-width buffer was created for modeling. Specifically, the order of the stream or river 

was multiplied by 100 m, and waterbodies larger than 1 km2 were buffered by 1 km. Thus, for modeling 

purposes, buffers were 400 m for fourth-order streams , 900 m for ninth-order rivers, and 1 km for 

waterbodies greater than 1 km2

. Third, all areas outside of buffer zones were masked so they would not 

be considered by the satellite model. These variable-width buffers adjacent to streams, rivers, and 

waterbodies represent the initial areas considered for modeling (fig. 4), comprising about 6 percent of 

the project area (59,887 km2

). 

Pinyon pine (Pinus edulis) and juniper (Juniperus spp.) were masked inside riparian buffer zones 

with high-resolution imagery because the satellite model commonly confuses them with riparian 

vegetation above 1,524-m elevation (Hatten and Sogge, 2007). Spectral confusion between pines and 

riparian vegetation usually is not a problem below 1,524-m elevation (Hatten and Paradzick, 2003), but 

for this project, all elevations were considered in order to include a small number of flycatcher sites in 

Colorado (Durst and others, 2008). Because high-elevation flycatcher habitat is problematic and results 

in more classification errors, the area of predicted flycatcher habitat was tracked above and below the 

1,524-m elevation boundary throughout the rangewide analysis, by management unit (fig. 4).

All urban and agricultural areas within proximity-to-water buffers were masked because

flycatcher territories do not occur in such areas (fig. 5). This processing operation involved several steps 

to accomplish. First, all urban and agricultural features identified in the Southwest and South-Central 

U.S. Geological Survey Gap Analysis Program (GAP) layers (Lowry and others, 2005) were extracted 

and masked. Second, the urban and agricultural mask was updated by digitizing omitted areas that were 

clearly visible on high-resolution digital aerial photographs. Omitted areas typically were agricultural 

fields that had increased in size in the last decade, or urban boundaries. Third, areas that were 

incorrectly classified as urban or agricultural were identified and the mask in such areas was removed, 

thereby including them in the modeling arena. 
9 

Figure 4. Map showing all fourth-order or larger streams and waterbodies, plus lakes greater than 1 square 

kilometer, and two elevation zones (1) below 1,524 m, and (2) above 1,524 m. 

The final masking step was an ocular assessment, whereby areas that were not organic in shape 

(for example, circles, lines, squares) were visually identified with the aid of high-resolution aerial orthophotographs and GIS, and excluded (fig. 6). This was primarily a refinement of the GAP layer because

agricultural areas oftentimes become spectrally confused with riparian areas along their interface. 

Although such areas were difficult to discern with 30-m resolution Landsat imagery, they were visible 

with 1-m resolution aerial photography, making their identification and digitization possible. 

Additionally, I excluded from the masks numerous agricultural fields that the Bureau of Reclamation 

(Reclamation) converted from agricultural plantations to native vegetation in hopes of providing 

flycatcher habitat; such areas retained a non-organic shape and were initially masked. Finally, all urban 

footprints that had enlarged since the last GAP assessment were identified and masked because they 

were not classified as urban at that time; such adjustments were necessary because flycatchers do not 

breed in urbanized areas.
10 

Figure 5. Map showing urban and agriculture areas identified and masked (that is, excluded) from the satellite 

model.
11 

Figure 6. Maps showing effects of masking on output from satellite model at confluence of the Colorado and Gila 

Rivers, Arizona and California. (A) No mask applied. (B) Distance-to-water mask applied. (C) Urban-agriculture 

mask applied. (D) Ocular mask applied. In these maps, the masks are additive from B to D.
12 

Application of the Satellite Model

Four steps were completed to project and test the satellite model rangewide (Hatten and 

Paradzick, 2003): 

1. Four GIS variables derived from Landsat Thematic Mapper (TM) imagery and 30-m resolution 

digital elevation models (DEMs) were created;

2. A logistical regression model was populated with four GIS variables and a probability grid was 

created;

3. Riparian vegetation was divided into probability classes; 

4. Model accuracy (sensitivity and specificity) was assessed within five probability classes and at a 

40-percent probability threshold, with flycatcher territory locations collected throughout 

flycatcher range. 

The satellite model uses a logistical regression equation developed by Hatten and Paradzick (2003) to 

calculate the probability of flycatcher habitat with the following equation: 

Habitat probability = exp(logit) / 1 + exp(logit) (1) 

where 

logit is 1.483(NDVI)+0.098(NDVIBEST)+0.034(FLOODPL)+0.648(NDVISTD)– 6.074, 

whereby the four variables are defined as follows:

1. NDVI = dense vegetation (NDVI>0.33) within a 30 × 30-m cell (0.09 ha); 

2. NDVIBEST = amount (percent) of densest vegetation (NDVI>0.41) within a 120-m radius 

(4.5-ha neighborhood); 

3. FLOODPL = amount (percent) of floodplain or flat terrain (<2.5 degrees) within a 360-m 

radius (41-ha neighborhood); 

4. NDVISTD = the standard deviation in NDVI (12 classes) within a 120-m radius (4.5-ha 

neighborhood). 

The three vegetation variables are extracted from Landsat TM imagery, and the FLOODPLAIN variable 

is extracted from a 30-meter resolution DEM. 

A set of custom GIS scripts were used to extract and process the information necessary to 

populate the satellite model (Hatten and Paradzick, 2003). Given the large size of the project area, each

Landsat scene was individually processed. Once completed, all 57 scenes were mosaicked together to 

form a continuous probability grid for classification and mapping purposes. The immense size of the 

project area required the division of the floodplain variable into 57 tiles that matched the Landsat 

footprints. As each Landsat scene was entered into the model, the script supplied the correct floodplain 

tile to complete the modeling process. Custom scripts also were used to reclassify the continuous 

probability grids into binary and five-probability-class format. Finally, the individual map tiles were 

mosaicked together to form a continuous habitat map for the entire range. The probability threshold was 

set at 40 percent to match the original habitat map (Hatten and Paradzick, 2003), but this threshold can 

easily be changed should the need arise. For example, if flycatchers occupied areas that were in a higher 

probability class (for example, 70 percent), increasing the threshold would reduce commission errors 

because less of the riparian zone would be considered suitable (Paxton and others, 2007). 
13 

Accuracy Assessment 

Habitat maps were created rangewide during 2013–15 and the fit and accuracy of the satellite 

model were assessed in several ways. Model fit was assessed by overlaying 758 flycatcher territory 

locations obtained in the summer of 2014 on five probability classes output by the satellite model with 

2014 Landsat imagery (fig. 7). The density of flycatcher territories was obtained by dividing territory 

numbers in each probability class by class area (Hatten and Paradzick, 2003). Past applications of the 

model have shown that the model is working correctly when flycatcher territory densities increase 

exponentially in high-probability classes. Model accuracy was assessed by calculating the percentage of 

flycatcher territories that occurred in predicted habitat at a 40-percent probability threshold because this 

was the optimum threshold determined in the original model. However, this threshold can be changed 

whenever model accuracy can be improved without increasing commission errors. 

Figure 7. Map showing locations of 758 Southwestern Willow Flycatcher (SWFL) territories that were used for 

model verification overlaid on 35 U.S. Fish and Wildlife Service management unit boundaries, Southwestern United 

States, summer 2014. Blue areas are all areas that the satellite model considered when predicting the occurrence 

of flycatcher habitat; white areas were excluded from modeling.
14 

Summarization and Display of Model Results

The broad scale of the project area and the large quantity of data output by the satellite model 

created challenges when characterizing and displaying information. The minimum mapping unit was an 

individual Landsat cell (0.09 ha), but there were more than 1 billion cells, so these challenges were 

handled hierarchically. Specifically, the habitat predictions were summarized at five different scales 

(mesh sizes). At the coarsest scale of analysis, the entire project area was considered irrespective of any 

internal management units. The remaining four analysis scales, in descending order of resolution, were 

6 State boundaries, 35 FWS management units, 88 FWS critical-habitat reaches, and 6,521 7.5-minute 

quadrangles. These five scales of resolution enabled a thorough analysis and display of predicted 

flycatcher habitat to satisfy many uses. For example, natural resource managers can review the areas of 

predicted flycatcher habitat within individual critical-habitat reaches to determine whether additional 

field visits or analysis are justified, or can determine if additional field surveys should be planned in 

locations that have seldom or ever been visited based on habitat maps developed from the 7.5-minute 

quadrangles. 

Habitat Change Detection

After a thorough GIS analysis of predicted flycatcher habitat was done during 2013–15, the 

results were presented in tabular and (or) graphical format, depending on the variability in the data and 

sample size. At the regionwide, State, and FWS management-unit scales, the total areas of predicted 

flycatcher habitat were computed and these results were shown in bar graphs and tables, allowing one to 

quickly compare which areas increased, decreased, remained the same, or fluctuated. The area of 

predicted flycatcher habitat that occurred inside 87 critical-habitat reaches was shown in tabular format 

because the variability in habitat and sample size made a bar graph impractical. In contrast, local and 

rangewide habitat change-detection maps were created, with results obtained from the 7.5-minutequadrangle analysis, as well as tabular output, because the 1:24,000 scale was a practical scale to work 

with. 

Tamarisk Leaf Beetle Impact Assessment

Three different modeling scenarios related to beetles and their associated effects on flycatcher 

habitat were examined. 

1. How beetles affected predicted flycatcher habitat along the beetle-infested lower Virgin 

River (southeast Nevada and northwestern Arizona), 

2. How beetles may affect predicted flycatcher habitat along the Lower Colorado River 

(southern Nevada and western Arizona), and, 

3. How beetles may affect predicted flycatcher habitat along the upper Gila River in 

southeastern Arizona. 

The Virgin River analysis differed from the beetle-impact simulations because it provided actual

information on how beetles affected flycatcher habitat from 2008 to 2015. The Virgin River case study 

also provided useful information when interpreting the two beetle-impact simulations done along the 

lower Colorado and upper Gila Rivers. 
15 

The two locations where beetle impacts on flycatcher habitat were simulated were based on the 

availability of high-quality vegetation maps that clearly delineated tamarisk patches (Bangle and others, 

2014; Neale and others, 2014). Tamarisk maps are necessary to do a beetle-impact simulation because

the beetles are tamarisk obligates that feed on the foliage of tamarisk, causing defoliation of the tree. 

Repeated defoliation of individual tamarisk trees can result in a severe dieback during the next season 

and ultimately mortality of the tree within several years (Deloach and Carruthers, 2004), affecting 

potential flycatcher breeding habitat. Beetle effects on predicted flycatcher habitat were simulated along 

the lower Colorado and upper Gila rivers in several steps. First, the vegetation maps were resampled to 

match the 30-m pixel resolution of the flycatcher maps (grids) with a 30-m majority filter (that is, the 

most dominant plant was assigned to a 30-m pixel). The cell centroids were snapped so that the pixels in 

the plant maps aligned with the flycatcher maps. Second, a conditional statement (CON) was used in 

ArcGISTM to reduce all Landsat NDVI values at tamarisk locations (cells) by 50 percent (value 

determined from Virgin River beetle analysis). Third, the flycatcher model was re-run using the altered 

NDVI values, which affected all three vegetation variables in the satellite model (see eq. 1). Fourth, a 

comparison was made on a cell-by-cell basis as to how predicted flycatcher habitat changed after the 

beetle-impact simulation with the ArcGISTM’ COMBINE function. Fifth, changes in predicted habitat 

were summarized across the entire reach and within distinct zones that had unambiguous start 

(upstream) and end (downstream) points. When comparing habitat maps created in different years, there 

were only four possible outcomes (Hatten and Paradzick, 2003): 

1. Habitat persisted (remained), 

2. Habitat disappeared (lost), 

3. Habitat formed (gained), or 

4. Unsuitable areas remained unsuitable. 

Habitat Time Series

An important aspect of the beetle analysis was to compare the magnitude of changes in 

flycatcher habitat to changes caused by natural stressors, such as drought, fire, or floods, as they all 

affect the quantity and quality of flycatcher habitat annually (Hatten and others, 2010; Orr and others, 

2014). Specifically, do beetle impacts on flycatcher habitat exceed the natural annual variability? If so, 

then by how much? If not, then perhaps the impacts to flycatchers are acceptable and can be considered 

as just one more natural stressor. Answering these questions required an historical record of flycatcher 

habitat for a given reach. To do so, two habitat time series were created, one for the lower Virgin River 

where beetle impacts have actually occurred, and one for the upper Gila River, where beetles have not 

arrived and impacts were simulated. In the lower Colorado River where beetles have not arrived, beetle 

effects on flycatcher habitat for 2015 were simulated without the aid of a habitat time series, with 

reliance instead on the knowledge gained from the Virgin River and upper Gila River habitat time series 

when interpreting results. 
16 

The two habitat time series were created with Landsat imagery dating from 1986 to 2015, 

obtained from the EROS satellite archive. Creating a habitat time series involved running the satellite 

model annually as close to the same date as possible and producing a habitat map for each respective 

year (Hatten and others, 2010). Habitat maps were created by applying a 40-percent probability 

threshold to the continuous probability grids output by the satellite model using the same methodology 

as the regionwide modeling effort. Finally, 95-percent confidence intervals (CIs) were calculated for the 

habitat time series and actual (Virgin River) and simulated (upper Gila River) habitat quantities were 

compared to the historical baseline. This approach allowed me to determine just how much of an impact 

beetles have had (lower Virgin River) and may have (upper Gila River) compared to the baseline 

variability caused by natural (for example drought, floods) and anthropogenic (for example, fire)

stressors. Given the size of the three modeled reaches (65–549 ha), they were divided into distinct zones 

in order to compare changes at a finer scale. 

Virgin River Beetle-Impact Analysis

The Virgin River beetle-impact analysis had three distinct objectives: 

1. To examine how beetles affected predicted flycatcher habitat upstream of Lake Mead where 

they have occurred since 2008 (Bateman and others, 2013); 

2. To determine how much NDVI values decreased in this reach so the same reduction factor 

could be applied to unaffected reaches in other parts of the State, thereby simulating beetle 

impacts for management purposes; and

3. To determine how much beetles affected flycatcher habitat by comparing the post-beetle 

habitat to a habitat time series from 1986 to 2015. 

For the first objective, a 65-km reach along the lower Virgin River was examined, from the 

Narrows (about 8 km upstream of Beaver Dam Wash) downstream to the Lake Mead National 

Recreation Area boundary. Changes in flycatcher habitat were assessed within seven zones separated by 

tributary junctions (fig. 8). For the second objective, habitat change detection (cell-by-cell) was done 

from 2010 to 2015 in an area severely impacted by beetles (Nagler and others, 2014). Although beetles

were first observed along the lower Virgin River in 2008, the habitat time series did not show beetle 

effects until 2011, which is consistent with previous research along this reach (Dobbs and others, 2012; 

Bateman and others, 2013; Bateman and Johnson, 2015), with radical habitat decreases occurring from

2013 to 2015. Thus, for change detection, 2010 was selected as the pre-beetle flycatcher habitat map 

and 2015 was selected as the post-beetle habitat map. For the second objective, I quantified how the 

mean NDVI values changed from 2010 to 2015 in areas the satellite model predicted suitable in 2010 

but were unsuitable in 2015, which I call an NDVI reduction factor. For the third objective, a habitat 

time series was created from 1986 to 2015, which required that two Landsat scenes be mosaicked 

together for each year since the study reach cut across two Landsat footprints (path rows 3934 and 

3935). The Landsat metadata associated with the Virgin River habitat time series are presented in 

appendix A, table A4. 
17 

Figure 8. Map showing boundary (polygon with seven segments) of lower Virgin River where a habitat time series 

was calculated and a change detection was done, Nevada and Arizona, 1986–2015. Seven zones were used to 

calculate changes along length of Virgin River.
18 

Lower Colorado River Beetle-Impact Simulation

The potential effects of beetles on predicted flycatcher habitat along 549 river kilometers (341 

river miles) of the Lower Colorado River were simulated, from US Route 93(0.5 km downstream of

Hoover Dam) to the Mexico border (fig. 9). This river section cuts across four Landsat scenes (pathrows 3935, 3936, 3836, 3837), so they were mosaicked together prior to running the satellite model 

(acquired in June and July 2015). This analysis differed from the Virgin and Gila River analyses in that 

the model was only run twice—once with unmodified imagery that represents pre-beetle riparian 

conditions, and once after the imagery was modified (that is, the NDVI reduction factor was applied 

before rerunning the model) to simulate beetle impacts on riparian vegetation. The FWS management 

unit boundaries were used to divide the lower Colorado River into four zones: Lower Gila, ParkerMexico, Bill Williams, and Hoover-Parker. 

Tamarisk locations were obtained for the beetle-impact simulation from a previous mapping 

effort that used 1-m resolution 2010 National Agriculture Imagery Program (NAIP) imagery, plus 

higher resolution imagery from 2011 (ESRI World Imagery) to assist with structural labels (Bangle and 

others, 2014). The authors used Trimble eCognition® to develop vegetation community maps similar to 

those developed along the Lower Colorado River in 1981 (Anderson and Ohmart, 1984). For 

consistency, Bangle and others (2014) selected minimum mapping units of 10 acres for 

cottonwood/willow/mesquite, 1 acre for marsh and water-related classes, and 25 acres for all other

vegetation types (for example, tamarisk). A lack of funding prevented an accuracy assessment of the 

vegetation map, but it is considered to be the highest-quality vegetation map to date for the study area. 
19 

Figure 9. Map showing Lower Colorado River boundary, tamarisk patches and structure (I–VI), U.S. Fish and 

Wildlife Service (FWS) Southwestern Willow Flycatcher (SWFL) management unit boundaries overlaid on a 

normalized difference vegetation index (NDVI) image obtained July 2015, California, Nevada, and Arizona.
20 

Upper Gila River Beetle-Impact Simulation

The potential effects of beetles on predicted flycatcher habitat were simulated along a 71-km 

reach of the upper Gila River, from the confluence of Bonita Creek downstream to Goodwin Wash 

(about 7 km downstream of Fort Thomas). The magnitude of change in predicted flycatcher habitat was 

examined after the beetle-impact simulation (using same method as described previously) by comparing 

results to the historical average and annual variability in flycatcher habitat. Specifically, a 30-year 

habitat time series was created on the upper Gila River (fig. 10), from 1986 to 2015, with Landsat 

imagery and the satellite model. Twenty-four of 29 Landsat scenes were obtained in June, with the 

remainder occurring in May, July, or August. The Landsat metadata associated with the upper Gila 

River habitat time series are presented in appendix A, table A5. Image quality was high for all dates, 

with the exception of 2012, which was omitted because of a lapse in Landsat 5 and Landsat 8 coverage. 

The natural variability in flycatcher habitat was calculated annually with a habitat time series, a bar 

graph was created of the time series, and a 95-percent CI was calculated for the historical record (1985–

2015). Finally, the post-beetle flycatcher habitat simulation was plotted and this result was compared to 

the historical habitat time series. 

The tamarisk locations used in the beetle-impact simulation are derived from an earlier 

vegetation mapping effort during late summer 2012 (Neale and others, 2014). High-resolution 

multispectral imagery (0.16 m) and 1-m resolution digital elevation data (obtained from lidar) were 

obtained on October 2, ortho-rectified, and co-registered. Field work was done in May 2013 to obtain 

comprehensive ground truth data, with dozens of locations visited. Trimble eCognition®, object‐based 

image analysis techniques were used to classify the remotely sensed data, sorting vegetation into 11 

classes:

1. Bare ground/sparse vegetation, 

2. Cottonwood/Gooddings willow, 

3. Defoliated vegetation, 

4. Grass/forb/low vegetation, 

5. Riparian grasses,

6. Riparian willow,

7. Shadow, 

8. Sparse shrub vegetation, 

9. Tamarisk,

10. Water, and

11. Willow.

Like the vegetation map for the Colorado River, the Gila River map was not assessed for accuracy, but 

it has been used extensively by the Gila Watershed Partnership for restoration activities and is 

considered to be the best map available (Orr and others, 2014). 
21 

Figure 10. Map showing project boundary, nine river reaches used in tamarisk leaf beetle-impact simulation, and 

nine vegetation classes, upper Gila River, Arizona, 1986–2015.
22 

Results

Satellite-Based Habitat Predictions

Rangewide

Rangewide predicted flycatcher habitat was 476.5 km2 (40-percent threshold) in 2013, 502.5 

km2 in 2014, and 756.5 km2 in 2015. Thus, predicted habitat was 5.5 percent more in 2014 than in 

2013,and 58.8 percent more in 2015. Summarized by State, the area of predicted habitat in 2015 (below 

1,524-m elevation), was least in Colorado, followed in ascending order by Utah, Nevada, Texas, 

California, Arizona, and New Mexico (fig. 11). This rank order remained almost the same when all 

elevations were considered, except for Colorado, where the predicted habitat was only 1.2 km2

 below 

1,524-m elevation, but 243.4 km2 above 1,524-m elevation. 

Figure 11. Graph showing area of predicted flycatcher breeding habitat below 1,524-meter elevation and for all 

elevations by State. Number above each bar is area of predicted habitat in square kilometers.
23 

U.S. Fish and Wildlife Service Management Units

Summarized by 35 FWS management units, the first, second, and third most predicted habitats 

below 1,524-m elevation were the Middle Rio Grande, Upper Gila, and Middle Gila/San Pedro 

management units in 2013 and 2014 (table 1; fig. 12). In 2015, this pattern was the same except for a 

large increase in predicted habitat in the Pecos management unit (New Mexico and Texas). Three 

management units are at elevations above the 1,524-m contour (Upper Rio Grande [New Mexico and 

Colorado], San Luis Valley [New Mexico and Colorado], and Dolores [Colorado]), without any lowelevation habitat. In contrast, 10 low-elevation management units (below 1,524-m elevation) have no 

high-elevation habitat. Predicted flycatcher habitat increased in 22 of 33 management units (66.7 

percent) below 1,524-m elevation from 2013 to 2015 (four units were excluded because they were too 

high), with the remainder showing a decrease (fig. 13). In contrast, predicted flycatcher habitat 

decreased in 6 of 28 management units (21.4 percent) at elevations above the 1,524-m contour (9 units 

were excluded because they were too low), and 21 of 27 (77.7 percent) increased in predicted habitat. 
24 

Figure 12. Graphs showing area of predicted flycatcher breeding habitat in U.S. Fish and Wildlife Service 

management units below and above 1,524-meter elevation, Southwestern United States, 2013–15. 
25 

Table 1. Changes in predicted flycatcher habitat in U.S. Fish and Wildlife Service management units, Southwestern 

United States, 2013–15.

[2013_a, 2014_a, and 2015_a: Predicted habitat (hectares) below 1,524-meter (m) elevation for each given year. 2013_c, 

2014_c, and 2015_c: Predicted habitat (hectares) above 1,524-m elevation for each given year. Δ_below: Change (in 

percent) below 1,524-m elevation. Δ_above: Change (in percent) above 1,524-m elevation. Data obtained with a satellite 

model at a 40-percent probability threshold. NA, not applicable] 

Management unit

Habitat area (hectares) Change (percent)

2013_a 2014_a 2015_a 2013_c 2014_c 2015_c Δ_below Δ_above

Texas Rio Grande 190.35 390.06 2,313.90 0.81 3.51 15.48 11.16 18.11

Lower Rio Grande 306.99 429.84 1,108.62 163.98 137.79 365.94 2.61 1.23

Santa Cruz 1,054.98 1,077.75 1,877.22 0.00 0.00 0.00 0.78 NA

Middle Gila/San Pedro 3,963.51 4,284.81 4,971.24 26.10 22.32 2.52 0.25 -0.90

Upper Gila 4,769.91 5,405.04 5,836.68 432.54 230.67 452.07 0.22 0.05

San Francisco 354.87 240.39 278.37 269.91 190.35 282.51 -0.22 0.05

Lower Gila 874.80 618.57 1,234.80 0.00 0.00 0.00 0.41 NA

Roosevelt 522.36 727.11 1,229.40 1,528.02 945.72 1,363.14 1.35 -0.11

San Diego 3,566.61 2,541.24 2,903.67 0.00 0.00 0.00 -0.19 NA

Parker-Mexico 2,655.18 3,383.46 2,897.10 0.72 0.36 0.00 0.09 -1.00

Hasayampa/Agua Fria 1,012.95 1,971.54 1,983.87 0.00 0.00 0.00 0.96 NA

Salton 1,086.66 874.71 957.15 52.56 58.86 57.69 -0.12 0.10

Santa Ana 2,600.37 2,371.41 2,371.14 16.83 12.69 63.27 -0.09 2.76

Middle Rio Grande 4,452.21 7,456.95 13,999.95 893.34 1,345.59 2,896.29 2.14 2.24

Bill Williams 1,762.65 1,662.57 1,748.34 5.67 8.46 8.28 -0.01 0.46

Santa Clara 2,535.21 1,823.85 1,562.40 71.55 81.18 72.54 -0.38 0.01

Verde 2,461.59 625.86 1,710.36 871.74 232.29 772.92 -0.31 -0.11

Pecos 502.56 1,940.31 12,269.25 474.03 1,396.98 6,513.30 23.41 12.74

Hoover-Parker 790.56 1,355.85 1,232.46 0.00 0.00 0.00 0.56 0.00

Upper Rio Grande 0.00 0.00 0.00 2,736.81 3,019.41 5,257.71 NA 0.92

Little Colorado 69.84 72.63 249.21 676.89 514.71 1,420.02 2.57 1.10

Mojave 388.08 582.84 600.30 32.31 34.74 31.05 0.55 -0.04

Santa Ynez 2,805.66 2,757.51 2,127.15 0.00 0.00 0.00 -0.24 0.00

San Luis Valley 0.00 0.00 0.00 6,256.62 7,323.84 9,772.29 NA 0.56

Upper San Juan 634.59 650.97 988.11 8,654.49 6,457.05 10,815.12 0.56 0.25

Middle Colorado 269.82 400.14 784.62 185.94 180.72 327.87 1.91 0.76

Lower San Juan 847.98 679.50 1,450.98 385.74 270.90 1,448.91 0.71 2.76

Dolores 0.00 0.00 0.00 3,939.21 2,058.30 6,211.62 NA 0.58

Virgin 374.13 361.53 450.63 251.55 226.35 256.59 0.20 0.02

Powell 42.21 138.42 147.42 75.24 71.91 126.09 2.49 0.68

Kern 3,363.84 2,338.38 2,423.97 4.14 3.15 6.93 -0.28 0.67

Sevier 0.00 0.00 0.00 550.89 666.27 1,140.39 NA 1.07

Amargosa 338.04 247.14 360.36 0.00 0.00 0.00 0.07 NA

Pahranagat 1,098.54 996.93 1,329.30 580.68 471.96 442.35 0.21 -0.24

Owen's 1,777.50 1,701.09 2,101.05 264.33 277.20 154.98 0.18 -0.41
26 

Figure 13. Graph showing changes in predicted flycatcher breeding habitat in U.S. Fish and Wildlife Service 

management units below and above 1,524-meter elevation and for all elevations, 2013–15.
27 

7.5-Minute Quadrangles

Summarized by 7.5-minute quadrangles, there was no predicted flycatcher habitat at any 

elevation from 2013 to 2015 in 69.4 percent (4,422) of the 6,521 quadrangles in the project area, but one 

must consider that masking excluded habitat predictions on streams smaller than 4th order or 

waterbodies less than 1 km2. These numbers changed when locations below 1,524-m elevation were 

examined, with no predicted habitat in 79.9 percent (5,212) of the quadrangles. The area of predicted 

flycatcher habitat ranged from 0.1 to 1,398 ha in a checkerboard pattern (fig. 14). Of the 50 top-ranked 

quadrangles below 1,524-m elevation (rank), predicted flycatcher habitat was 42 percent (21) in 

Arizona, 32 percent (16) in New Mexico, 22 percent (11) in California, and 2 percent (1) in Utah and 

Nevada, respectively (fig. 15; table 2). The 10 top-ranked quadrangles in descending order were Paraje 

Well and San Marcial (New Mexico; fig. 16A), San Carlos Reservoir (southeastern Arizona; fig. 16B), 

Corona North (California), Romero Canyon (New Mexico), Lake Mcmillan North (New Mexico), Prado 

Dam (California), Weldon (California), Needles (Arizona), and Dewey Flat (Arizona) (table 2). 
28 

Figure 14. Map showing area of predicted flycatcher breeding habitat (averaged across 2013–15) at all elevations 

in 6,521 U.S. Geological Survey 7.5-minute quadrangles, Southwestern United States. Only areas within a 

prescribed distance of streams and greater than stream-order 3, or within 1 kilometer of a lake or reservoir, were 

included in this analysis. Areas outside the selection criteria are excluded (white), and areas with no habitat are 0 

(red).
29 

Figure 15. Graph showing changes in predicted flycatcher habitat for the 50 top-ranked U.S. Geological Survey 

7.5-minute quadrangles below 1,524-meter elevation, 2014–15. Predicted flycatcher habitat increased in 34 

quadrangles (68 percent) and decreased in 16 quadrangles (32 percent). 
30 

Figure 16. Maps showing predicted flycatcher breeding habitat at (A) Paraje Well and San Marcial, New Mexico, 

June 2013, and (B) San Carlos Reservoir, south-central Arizona. The 7.5-minute U.S. Geological Survey 

quadrangles and their names are overlaid.
31 

Table 2. Top-ranked 50 7.5-minute U.S. Geological Survey quadrangles by average area (hectares) of predicted 

flycatcher habitat, Southwestern United States, 2013–15.

[Habitat estimates were obtained from the satellite model at a 40-percent probability threshold. USGS_QD_ID: Quadrangle 

identification number. QUAD_NAME: Quadrangle name. ST_NAME State in which quadrangle occurs. 3yravg_p40:

Average area (hectares) of flycatcher habitat in quadrangles below 1,524-meter elevation shown in descending order]

USGS_QD_ID QUAD_NAME ST_NAME 3yravg_p40

33107-E1 Paraje Well New Mexico 1,398.27

33106-F8 San Marcial New Mexico 1,331.97

33110-B4 San Carlos Reservoir Arizona 1,030.47

33117-H5 Corona North California 943.53

33107-D2 Romero Canyon New Mexico 907.8

32104-F3 Lake McMillan North New Mexico 759.87

33117-H6 Prado Dam California 655.98

35118-F3 Weldon California 633.33

34114-G5 Needles Arizona 585.69

33110-B3 Dewey Flat Arizona 578.82

34113-C5 Artillery Peak Arizona 569.4

33107-C2 Black Bluffs New Mexico 569.31

33110-B2 Calva Arizona 567.03

33106-G7 San Antonio SE New Mexico 536.85

34106-H6 Isleta New Mexico 534.27

37115-D2 Ash Springs Nevada 513.57

33114-C6 Cibola Arizona 470.46

32110-F5 Clark Ranch Arizona 460.38

33109-A8 Fort Thomas Arizona 447.78

34114-C1 Monkeys Head Arizona 440.61

33112-C6 Hassayampa Arizona 427.89

34106-G6 Los Lunas New Mexico 416.19

34106-C7 La Joya New Mexico 404.94

33114-F5 Blythe NE California 400.02

37118-D4 Fish Slough California 376.47

31111-E1 Tubac Arizona 375.63

32109-H8 Eden Arizona 372.12

34119-D3 Matilija California 365.16

32110-H6 Dudleyville Arizona 363.15

33112-C5 Buckeye Arizona 358.2

33106-H7 San Antonio New Mexico 337.56

34111-A6 Chalk Mountain Arizona 320.43

33112-C2 Laveen Arizona 315.87

33117-C3 Morro Hill California 313.74

32104-G3 Spring Lake New Mexico 312.51

32110-G6 Lookout Mountain Arizona 308.22
32 

USGS_QD_ID QUAD_NAME ST_NAME 3yravg_p40

37110-B5 No Mans Mesa North Utah 306.63

34113-F5 Wikieup Arizona 304.53

34120-G5 Casmalia California 297.15

34106-E7 Veguita New Mexico 295.68

32110-E5 Peppersauce Wash Arizona 294.12

33107-F1 Fort Craig New Mexico 288.18

33114-D6 Palo Verde California 287.91

33114-B6 Picacho NW California 287.7

34106-F6 Tome New Mexico 284.49

36119-G4 Piedra California 276.42

32108-H5 Cliff New Mexico 272.91

33110-C4 San Carlos Arizona 272.52

34106-D7 Abeytas New Mexico 267.09

33110-A1 Geronimo Arizona 265.83

U.S. Fish and Wildlife Service Critical-Habitat Reaches

Predicted flycatcher habitat was the most in the Rio Grande recovery unit (table 3) when 

averaged across 2013–15 (8,854 ha), followed in descending order by the Gila (6,022 ha), Coastal 

California (2,911 ha), Little Colorado (1,416 ha), Basin and Mojave (1,258 ha), and Upper Colorado 

(258 ha). The top five critical-habitat reaches, ranked in descending order for predicted flycatcher

habitat (averaged across 2013–15) were Middle Rio Grande (7,165 ha), San Pedro River (2,003 ha), 

Gila River–mid San Carlos (1,247 ha), Rio Grande–San Luis Valley National Wildlife Refuge (1,107

ha), and Santa Clara River (962 ha). Four critical-habitat reaches contained no habitat in any year: Mill 

Creek–Salton, Paria River, Bautista Creek–east, and Bautista Creek–west.
33 

Table 3. Predicted flycatcher habitat in 87 U.S. Fish and Wildlife Service critical-habitat reaches (Unit name), 

located within 6 recovery units and 35 management units, Southwestern United States, 2013–15.

[Critical-habitat reach locations are shown in figure 1 by reach identification number (ID No.). Predicted habitat was 

calculated at a 40-percent probability threshold, as determined from the satellite model. Mill Creek has the same ID No. for 

two management units because it was split in two in the recovery plan]

Recovery unit Management unit Unit name

Habitat area (hectares) ID

2013 2014 2015 No.

Basin and Mojave Amargosa Amargosa River 11.6 2.5 12.2 25

Basin and Mojave Amargosa Ash Meadows NWR 0.2 0.2 0.2 56

Basin and Mojave Amargosa Willow Creek 10.5 3.0 10.4 44

Basin and Mojave Kern Canebrake Creek 10.1 8.8 9.5 43

Basin and Mojave Kern Kern River South Fork 891.2 646.0 559.8 13

Basin and Mojave Mojave Deep Creek 0.0 0.0 2.6 48

Basin and Mojave Mojave Holcomb Creek 5.4 3.9 6.8 47

Basin and Mojave Mojave Mojave River 357.4 499.1 533.6 11

Basin and Mojave Salton Mill Creek 0.0 0.0 0.0 79

Basin and Mojave Salton San Felipe Creek 77.4 51.9 60.3 80

Coastal California San Diego Agua Hedionda Creek - east 12.5 7.7 8.4 84

Coastal California San Diego De Luz Creek 0.3 0.0 0.0 86

Coastal California San Diego Pilgrim Creek 7.3 7.3 7.1 87

Coastal California San Diego San Diego River El Capitan - north 11.3 5.6 10.2 72

Coastal California San Diego San Diego River El Capitan - south 28.3 6.8 9.9 74

Coastal California San Diego San Luis Rey River - east 4 0.6 1.3 0.7 67

Coastal California San Diego San Luis Rey River - east 5 1.3 0.1 1.8 70

Coastal California San Diego San Luis Rey River - east 6 0.3 0.3 1.7 69

Coastal California San Diego San Luis Rey River - east 7 45.6 38.4 40.1 68

Coastal California San Diego San Luis Rey River - middle 3 1.5 2.1 1.2 65

Coastal California San Diego San Luis Rey River - west 1 378.9 256.1 265.1 66

Coastal California San Diego San Luis Rey River - west 2 43.8 36.1 30.7 88

Coastal California San Diego Santa Margarita River 82.8 71.2 87.1 64

Coastal California San Diego Santa Ysabel Creek - east 27.1 16.6 31.8 85

Coastal California San Diego Sweetwater River 53.3 48.7 47.1 75

Coastal California San Diego Temecula Creek 6.5 3.2 3.4 82

Coastal California San Diego Temescal Creek 13.1 6.6 7.9 71

Coastal California Santa Ana Bautista Creek - east 0.0 0.0 0.0 78

Coastal California Santa Ana Bautista Creek - middle 1.0 0.4 1.6 77

Coastal California Santa Ana Bautista Creek - west 0.0 0.0 0.0 76

Coastal California Santa Ana Bear Creek 40.8 42.0 48.2 60

Coastal California Santa Ana Mill Creek 29.3 30.0 32.8 79

Coastal California Santa Ana Oak Glen Creek 15.9 14.6 12.3 63

Coastal California Santa Ana San Timoteo Creek 18.4 17.6 23.1 81

Coastal California Santa Ana Santa Ana River - east 50.1 44.6 55.1 62

Coastal California Santa Ana Santa Ana River - middle 5.1 21.5 28.4 83

Coastal California Santa Ana Santa Ana River - west 528.3 510.3 499.8 59

Coastal California Santa Ana Waterman Creek 10.2 9.8 12.7 61

Coastal California Santa Clara Big Tujunga Canyon 158.1 130.1 136.1 26

Coastal California Santa Clara Castaic Creek 7.9 0.7 0.8 46

Coastal California Santa Clara Piru Creek 10.7 41.0 61.3 27

Coastal California Santa Clara Santa Clara River 1,191.9 832.8 860.6 9

Coastal California Santa Clara Ventura River 134.4 81.7 95.3 30

Coastal California Santa Ynez Mono Creek 25.2 19.0 1.7 45
34 

Recovery unit Management unit Unit name

Habitat area (hectares) ID

2013 2014 2015 No.

Coastal California Santa Ynez Santa Ynez River - east 97.4 89.1 60.9 29

Coastal California Santa Ynez Santa Ynez River - middle 0.6 0.4 0.0 28

Coastal California Santa Ynez Santa Ynez River - west 369.1 306.3 139.5 10

Gila Hassayampa and Agua Fria Hassayampa River 63.9 43.7 1.4 23

Gila Middle Gila and San Pedro Gila River - San Pedro 148.3 670.4 792.1 38

Gila Middle Gila and San Pedro San Pedro River 2,055.3 2,014.9 1,939.8 2

Gila Roosevelt Salt River 17.9 115.7 305.6 37

Gila Roosevelt Tonto Creek 120.2 145.6 233.6 36

Gila San Francisco San Francisco River - east 76.1 48.2 114.0 8 

Gila San Francisco San Francisco River - NE 2.3 10.0 2.5 1

Gila San Francisco San Francisco River - NW 7.3 23.9 0.0 49

Gila San Francisco San Francisco River - west 50.5 30.2 13.5 19

Gila Santa Cruz Cienega Creek 104.9 97.0 89.0 21

Gila Santa Cruz Santa Cruz River 580.5 621.1 777.9 22

Gila Upper Gila Gila River - east 169.1 258.3 313.6 32

Gila Upper Gila Gila River - mid San Carlos 1,069.5 1,239.2 1,432.1 33

Gila Upper Gila Gila River - south 275.0 386.9 395.5 3 

Gila Verde Verde River - middle N 46.6 50.9 56.2 54

Gila Verde Verde River - middle S 43.7 39.1 7.3 35

Gila Verde Verde River - upper M 216.5 7.7 85.1 5 

Gila Verde Verde River - upper N 288.5 2.2 151.7 4 

Gila Verde Verde River - upper S 115.0 15.3 85.1 6

Little Colorado Bill Williams Big Sandy River - Alamo 145.9 144.0 168.0 7 

Little Colorado Bill Williams Big Sandy River - north 345.2 334.3 328.1 40

Little Colorado Bill Williams Bill Williams River - Alamo 268.0 271.8 276.3 41

Little Colorado Bill Williams Bill Williams River - Lincoln 42.3 37.5 42.8 24

Little Colorado Bill Williams Bill Williams River - NWR upper 220.1 186.2 146.8 55

Little Colorado Bill Williams Santa Maria River 225.0 220.9 227.1 42

Little Colorado Little Colorado Little Colorado River 40.6 12.3 35.1 50

Little Colorado Little Colorado Little Colorado River West Fork 0.4 0.0 0.2 12

Little Colorado Pahranagat Pahranagat - NWR 59.4 46.1 67.8 39

Little Colorado Virgin Virgin River 119.9 109.6 126.7 34

Rio Grande Middle Rio Grande Rio Grande - middle 3,578.0 6,052.5 11,864.6 14

Rio Grande San Luis Valley Conejos River SLV 230.5 229.1 348.6 57

Rio Grande San Luis Valley Rio Grande SLV - NWR 1,353.1 703.0 1,266.2 15

Rio Grande San Luis Valley Rio Grande SLV - south 1.3 4.6 46.8 51

Rio Grande Upper Rio Grande Coyote Creek 34.5 44.1 63.9 16

Rio Grande Upper Rio Grande Rio Fernando de Taos 1.0 1.2 2.3 31

Rio Grande Upper Rio Grande Rio Grande - Between Pueblos 6.0 8.3 14.1 52

Rio Grande Upper Rio Grande Rio Grande - upper 192.6 194.7 260.5 18

Rio Grande Upper Rio Grande Rio Grande del Rancho 9.4 9.5 42.8 17

Upper Colorado Powell Paria River 0.0 0.0 0.0 53

Upper Colorado San Juan Los Pinos River 168.9 161.4 197.7 58

Upper Colorado San Juan San Juan River - west 75.2 14.8 157.2 20
35 

Model Accuracy

At a 40-percent cutpoint (binary habitat map), 88 percent (667) of flycatcher territories were

inside predicted habitat in 2014, whereas 12 percent (91) were outside predicted habitat (omission). Of 

the 91 territories outside predicted habitat (less than 40-percent probability), 90 percent of them (82) 

were in the first or second probability classes (1–39-percent probabilities), whereas 10 percent were

outside any predicted habitat (0 percent). The Euclidean function in ArcGISTM indicated that the average 

distance of a flycatcher territory from predicted habitat was 74 m (standard deviation of 74 m), and the 

maximum distance was 277 m (about 9 cells). The density of flycatcher territories in 2014 was 

exponentially related to the five model probability classes (fig. 17), with 0.0005 territories per hectare in 

class 1, and 0.0177 territories per hectare in class 5. More than 99 percent of the variability in territory 

densities among the five density classes is explained in equation 2, indicating an excellent fit between 

model probabilities and flycatcher territories: 

D = 0.0002e0.9114x (2) 

where 

 D is flycatcher territory density (number of territories per hectare), and

X is probability class. 

The distributions of flycatcher territories are presented at five locations, overlaid on binary (40 

percent probability threshold) and five-probability maps: 

1. Owens River (fig. 18A–18B),

2. Alamo Lake (fig. 19A–19B),

3. Rio Grande-Elephant Butte Reservoir (fig. 20A–20B),

4. Colorado River-Havasu National Wildlife Refuge (fig. 21A–21B), and 

5. Santa Margarita and San Luis Rey Rivers (fig. 22A–22B). 

Many other areas were equally compelling to visualize, but were omitted for space considerations.
36 

Figure 17. Graph showing density of flycatcher territories (n = 747) in 2014 in five probability classes output by the 

satellite model.
37 

Figure 18. Flycatcher habitat suitability maps of the upper Owens River, California. (A) Produced by the satellite 

model at a 40-percent probability threshold; cells with values greater than 40 percent are considered suitable and 

cells with values less than 40 percent are considered unsuitable. (B) Produced at 20-percent probability thresholds; 

probability class 1 is the lowest ranked habitat and class 5 is the highest. Southwestern Willow Flycatcher (SWFL) 

territories were located in summer 2014 by the California Department of Fish and Wildlife.
38 

Figure 18.—Continued
39 

Figure 19. Flycatcher habitat suitability maps of Alamo Lake, Arizona. (A) Produced by the satellite model at a 40-

percent probability threshold; cells with values greater than 40 percent are considered suitable and cells with 

values less than 40 percent are considered unsuitable. (B) Produced at 20-percent probability thresholds; 

probability class 1 is the lowest ranked habitat and class 5 is the highest. Southwestern Willow Flycatcher (SWFL) 

territories were located in summer 2014 by the Bureau of Reclamation. 
40 

Figure 19.—Continued
41 

Figure 20. Flycatcher habitat suitability maps of Rio Grande (Elephant Butte Reservoir), New Mexico. (A) 

Produced by the satellite model at a 40-percent probability threshold; cells with values greater than 40 percent are 

considered suitable and cells with values less than 40 percent are considered unsuitable. (B) Produced at 20-

percent probability thresholds; probability class 1 is the lowest ranked habitat and class 5 is the highest.

Southwestern Willow Flycatcher (SWFL) territories were located in summer 2014 by the Bureau of Reclamation.
42 

Figure 20.—Continued
43 

Figure 21. Fflycatcher habitat suitability maps of Havasu National Wildlife Refuge (Topock Marsh), Arizona and 

California. (A) Produced by the satellite model at a 40-percent probability threshold; cells with values greater than 

40 percent are considered suitable and cells with values less than 40 percent are considered unsuitable. (B) 

Produced at 20-percent probability thresholds; probability class 1 is the lowest ranked habitat and class 5 is the 

highest. Southwestern Willow Flycatcher (SWFL) territories were located in summer 2014 by the Bureau of 

Reclamation.
44 

Figure 21.—Continued
45 

Figure 22. Flycatcher habitat suitability maps of the Santa Margarita and San Luis Rey Rivers (Camp Pendleton), 

California. (A) Produced by the satellite model at a 40-percent probability threshold; cells with values greater than 

40 percent are considered suitable and cells with values less than 40 percent are considered unsuitable. (B) 

Produced at 20-percent probability thresholds; probability class 1 is the lowest ranked habitat and class 5 is the 

highest. Southwestern Willow Flycatcher (SWFL) territories were located in summer 2014 by the U.S. Geological 

Survey. 
46 

Figure 22.—Continued
47 

Omission errors were acceptable in all areas examined, without any identified performance 

issues by geographic area, but commission errors were high wherever the model was applied. 

Commission errors are locations (30 × 30-m cells) determined through the model to be suitable but that 

did not contain a flycatcher territory (Story and Congalton, 1986). Unlike omission errors, which really 

are classification errors, commission errors do not provide proof that the model is wrong because a 

location may actually be suitable but not contain a flycatcher territory. Given the scarcity of flycatchers 

throughout their range, it is not surprising that many locations considered suitable by the model were 

empty. The reader may find a more detailed discussion related to model accuracy in section, 

“Discussion.” 

Habitat Change Detection

Predicted flycatcher habitat changed directionally from 2013 to 2015, with habitat generally 

increasing eastwards and decreasing westwards (fig. 23). Specifically, predicted habitat increased 

substantially in New Mexico and Texas; Arizona, Utah, and Nevada were mixed (both increases and 

decreases); and California was largely decreases. Habitat increases outnumbered decreases, with 68 

percent of the top-ranking quadrangles increased and 32 percent decreased. The top-ranking 

quadrangles were in Texas (48 percent), followed by New Mexico (36 percent), Arizona (14 percent), 

and California (2 percent). In contrast, the 50 quadrangles with the largest decreases in predicted habitat 

were in California (56 percent), Arizona (32 percent), Nevada (8 percent), and New Mexico (4 percent). 

Change in predicted flycatcher habitat was greatest at Lake Mcmillan, Texas, from June 2014 to June 

2015, increasing by 132 percent in the Lake Mcmillan North quadrangle, and 143 percent in the Spring 

Lake quadrangle (fig. 24). The second largest changes were observed near Fort Craig and San Marcial, 

New Mexico, where predicted habitat increased 177 percent in the Fort Craig quadrangle and 98 percent

in the San Marcial quadrangle (fig. 25).
48 

Figure 23. Map showing relative changes in predicted flycatcher habitat, as determined from a satellite model at a 

40-percent probability threshold, Southwestern United States, 2013–15.
49 

Figure 24. Map showing changes in predicted flycatcher habitat at Lake McMillan, Texas, June 2014–June 2015. 

Predicted habitat in 2014 completely overlays 2015 predicted habitat, which produces a hybrid color (brown). 

Predicted habitat increased 132 percent in the Lake McMillan North quadrangle, and 143 percent in the Spring 

Lake quadrangle.
50 

Figure 25. Map showing changes in predicted flycatcher habitat in Fort Craig and San Marcial quadrangles, New 

Mexico, 2014–15. Predicted habitat in 2014 completely overlays 2015 predicted habitat, which produces a hybrid 

color (brown). Predicted habitat increased 177 percent in the Fort Craig quadrangle and 98 percent in the San 

Marcial quadrangle between years.
51 

Tamarisk Leaf Beetle Impact Assessment

Lower Virgin River

The Virgin River habitat time series produced a cyclical pattern in predicted flycatcher habitat

from 1986 to 2015 (fig. 26), with a yearly average of 700.3 ha and a standard deviation of 175.6 ha (95-

percent CI [356 ha, 1,044 ha]). When calculating the historical average, I did not include habitat for 

2013–15 because it was well outside the 95-percent CIs observed from 1986 to 2012. Given this large 

deviation, which was clearly the result of beetles, I compared the last years (2013–15) to the historical 

average (1986 and 2012). Although the overall pattern from 1986 to 2015 was variable, a clear pattern 

of habitat formation, persistence, and decreases were evident. For descriptive purposes, the time series 

was divided into five distinct periods characterized by low, moderate, and high habitat production. 

Specifically, the period from 1986 to 1994 produced moderate amounts of habitat (503 and 710 ha), 

with the exception of 1992, which only produced 355 ha. The period from 1995 to 2001 was the most 

robust period of habitat production, when 3 of the top 5 years exceeded 900 ha. The period from 2002 to

2007 was another moderate period, producing from 540 to 749 ha per year. The years 2008–2010 were

another robust period, producing from 867 to 987 ha per year, with 2 of the top 5 years (912 and 937 

ha). The year 2011 was a transition year that separated a robust period (2008–2010) from an 

exceptionally poor period. Specifically, habitat was less than 100 ha per year (Ῡ = 88 ha) for 2013–15, 

falling well short of the lower 95-percent CI (356 ha), with just 12.6 percent of the historical average (-

3.5 standard deviations). 

Comparing the years 2010 and 2015, which are pre- and post-beetle years, respectively, there 

were 937 ha of predicted habitat in 2010, and 87 ha in 2015, representing a 90.7 percent loss (fig. 27). 

Decomposed further, 94.19 percent of predicted habitat was lost, 5.81 percent persisted, and there was a 

3.53 percent gain in new habitat (fig. 28). The NDVI values decreased by 48.9 percent where predicted 

flycatcher habitat was lost, and by 9.5 percent where habitat persisted (fig. 29). In areas that were 

predicted to be unsuitable both in 2010 and 2015, NDVI values actually increased by 6.7 percent. 

Finally, in regions that became suitable in 2015, NDVI values increased by 101.1 percent. Dividing the 

channel into seven change zones (fig. 30), zones 1–6 lost from 28.7 to 100.0 percent of predicted 

habitat, whereas zone 7 showed an increase of 160.0 percent (table 4). Predicted habitat was not evenly 

distributed in these zones, making the impacts disproportionate. For example, zone 1 (Lake Mead 

National Recreation Area boundary to Halfway Wash) had 618.3 ha of predicted habitat in 2010, which 

is over five times more than any other zone, and it lost 99.6 percent of its habitat. In contrast, zone 7 

(Beaver Dam Wash to the Narrows) had only 0.5 ha of predicted habitat in 2010 and 1.2 ha of predicted 

habitat in 2015, marking a 140 percent increase, but it only amounted to 0.8 ha. 
52 

Table 4. Area of predicted flycatcher habitat before and after tamarisk leaf beetle-impact analysis was done for the 

lower Virgin River, Nevada and Arizona, 2010 and 2015. 

[The satellite model produced the habitat estimates at a 40-percent probability threshold. Zone: There are seven lower Virgin 

River zones, as shown in figure 27. Habitat_pre: Flycatcher habitat before leaf beetle-impact analysis was done.

Habitat_post: Flycatcher habitat after leaf beetle-impact analysis was done. Change: Changes in predicted flycatcher 

habitat from 2010 (pre-beetle) to 2015 (post-beetle). NHDPlus tributary junctions: The seven zones are delineated by 

tributary junctions available in the U.S. Geological Survey NHDPlus dataset]

Zone

Habitat area (hectares) Change

(percent) NHDPlus tributary junctions Habitat_pre Habitat_post

1 618.3 2.6 -99.6 Lake Mead National Recreation Area boundary to 

Halfway Wash

2 51.9 8.4 -83.9 Halfway Wash to Nickel Creek

3 13.1 0.0 -100.0 Nickel Creek to Toquop Wash

4 100.7 42.8 -57.6 Toquop Wash to Pulsipher Wash

5 116.5 6.4 -94.5 Pulsipher Wash to Sand Hollow Wash

6 36.7 26.2 -28.7 Sand Hollow Wash to Beaver Dam Wash

7 0.5 1.2 160.0 Beaver Dam Wash to the Narrows (about 8 kilometers 

upstream of Beaver Dam Wash)
53 

Figure 26. A habitat time series of the lower Virgin River (see fig. 8), Nevada and Arizona, produced by running 

the satellite model repeatedly from 1986 to 2015 (2012 omitted). Mean value for 2013–15 was 88 ha, or 12.6 

percent of historical average (1986–2011). See table A4 for image acquisition dates for each year. 
54 

Figure 27. Map showing changes in predicted flycatcher habitat along the lower Virgin River, Nevada and Arizona, 

2010–2015, as determined from a satellite model at a 40-percent probability threshold.
55 

Figure 28. Graph showing changes in area of predicted flycatcher habitat in the lower Virgin River, Nevada and 

Arizona (fig. 27), 2010–2015. Specifically, 94.19 percent of habitat was lost, 5.81 percent persisted, and 3.53 

percent of new habitat was gained. Satellite model was used to predict flycatcher habitat at a 40-percent probability 

threshold.

Figure 29. Graph showing changes in NDVI values in four zones from 2010 (June 22) to 2015 (June 25—) in the 

lower Virgin River, Nevada and Arizona. Unsuitable (was not suitable flycatcher habitat in 2010 or 2015); Declined 

(flycatcher habitat was lost); Persisted (habitat persisted); and Gained (habitat formed). Thus, NDVI values 

declined 48.9 percent in areas where flycatcher habitat was lost (Declined) from 2010 to 2015, while NDVI 

increased 101.1 percent in locations where flycatcher habitat formed.
56 

Figure 30. Graph showing changes in predicted flycatcher habitat inside seven zones (fig. 27) along the lower 

Virgin River, Nevada and Arizona, 2010–2015. A satellite model was used to predict flycatcher habitat at a 40-

percent probability threshold.

Lower Colorado River

A beetle-impact simulation for the lower 549 km of the Colorado River resulted in a 33.5-

percent decrease in predicted flycatcher habitat (table 5), with 4,658 ha of predicted flycatcher habitat 

before the simulation and 3,099 ha afterwards. The Lower Colorado River reach was too large to show 

simulated changes efficiently, so close ups were provided of Topock Marsh (fig. 31) and Cibola 

National Wildlife Refuge (fig. 32). A close examination of flycatcher territories at Topock Marsh in 

2014 indicated that 79 percent of the territories occurred in habitat that is predicted to persist should 

beetles arrive, whereas 10.3 percent occurred in habitat predicted to be lost. The remaining three 

territories were not found in predicted habitat and are thus classified as omission error. The four 

Reclamation conservation areas did not have much predicted flycatcher habitat in 2015, so the beetleimpact simulation also resulted in small changes in these areas. In contrast, the Cibola National Wildlife 

Refuge showed much more promising results in the Reclamation conservation zones, with a substantial 

amount of predicted flycatcher habitat persisting after the beetle-impact simulation in numerous units. 
57 

Table 5. Changes in predicted Southwestern Willow Flycatcher habitat along lower Colorado River, California, 

Nevada, and Arizona, after a beetle-impact simulation was done in 2015.

[The satellite model produced the habitat estimates at a 40-percent probability threshold. Management Unit: Amount of 

predicted flycatcher habitat is divided into U.S. Fish and Wildlife Service flycatcher management units (see fig. 9). 

Pre_beetle is predicted flycatcher habitat (hectares) before a leaf beetle-impact simulation was done; Post_beetle is 

predicted flycatcher habitat (hectares) after a leaf beetle-impact simulation was done] 

Management unit

Habitat area (hectares) Change 

Pre_beetle Post_beetle (percent)

Lower Gila 69 37 -46.7

Parker-Mexico 2,899 1,843 -36.4

Bill Williams 557 510 -8.5

Hoover-Parker 1,133 710 -37.3

Total 4,658 3,099 -33.5

Divided into three distinct habitat states, the beetle-impact simulation along the Lower Colorado 

River resulted in 2.4 percent new habitat (113.5 ha), 35.9 percent lost habitat (1,672 ha), and 64.1 

percent (2,986 ha) habitat that persisted (fig. 33). The Lower Gila FWS flycatcher management unit 

(fig. 9) is predicted to lose 46.7 percent habitat, Parker-Mexico to lose 36.4 percent, Bill Williams to 

lose 8.5 percent, and Hoover-Parker to lose 37.3 percent (fig. 34). The magnitude of loss indicates the 

proportion of tamarisk that comprised (overlapped) predicted flycatcher habitat, with the Bill Williams 

having the least amount of overlap and the Lower Gila having the most amount of overlap. 
58 

Figure 31. Map showing simulation results showing how tamarisk leaf beetles may affect flycatcher habitat at 

northern part (Topock Marsh) of Havasu National Wildlife Refuge along lower Colorado River, Arizona and 

California. Three possible change outcomes are (1) habitat persisted, (2) habitat declined (lost), and (3) habitat 

gained. Bureau of Reclamation conservation area boundaries are overlaid. Satellite model was used to predict 

flycatcher habitat at a 40-percent probability threshold.
59 

Figure 32. Map showing simulation results showing how tamarisk leaf beetles may affect flycatcher habitat at 

Cibola National Wildlife Refuge along lower Colorado River, Arizona and California. Three possible change 

outcomes are (1) habitat persisted, (2) habitat declined (lost), and (3) habitat formed (gained). Bureau of 

Reclamation conservation boundary units are overlaid. Flycatcher model was used to predict flycatcher habitat at a 

40-percent probability threshold.
60 

Figure 33. Graph showing three states of flycatcher habitat related to tamarisk leaf beetle-impact simulation along 

lower Colorado River, Arizona and California (see fig. 9). There were 4,658 ha of predicted habitat before the 

simulation, and 3,099 ha post-beetle simulation, representing a 33.5 percent reduction. Specifically, the postsimulation produced 113.5 ha of new flycatcher habitat, while 2,985.7 ha persisted and 1,672.3 ha were lost. A 

satellite model was used to predict flycatcher habitat at a 40-percent probability threshold.

Figure 34. Graph showing changes in predicted flycatcher habitat along the lower Colorado River, by U.S. Fish 

and Wildlife Service management units(see fig. 9), after simulating a tamarisk leaf beetle infestation, Arizona and 

California. This analysis required that I run the flycatcher satellite model twice: once before the arrival of tamarisk 

leaf beetles (using 2015 imagery), and secondly, after NDVI values were artificially lowered in areas where 

tamarisk currently occurs. A satellite model was used to predict flycatcher habitat at a 40-percent probability 

threshold.
61 

Upper Gila River

The upper Gila River habitat time series indicated a cyclical pattern in predicted flycatcher 

habitat from1986 to 2015 (fig. 35), with a mean of 1,214 ha per year, a standard deviation of 277 ha, 

and 95-percent CIs of 672–1,756 ha. The year 2002 produced the least amount of predicted habitat, at 

684 ha, whereas 2008 produced the greatest amount (1,850 ha). Although the time series is ragged in 

appearance, there was a clear bimodal shape—1986–98 comprised the first hump, producing from 880 

to 1,746 ha; 2007–15 comprised the second hump, producing from 936 to 1,850 ha; and 1999–06 

comprised a trough, producing from 684 to 1,323 ha. The only year that exceeded the 95-percent CIs

was 2008 (1,850 ha), but 1990 almost exceeded it (1,746 ha), whereas 2002 almost fell below it (684 

ha).

A beetle-impact simulation for a 71-km reach of the upper Gila River resulted in a 53.1-percent

loss in predicted flycatcher habitat, with 1,308 ha before the simulation and 613 ha afterwards, 

decreasing below two standard deviations calculated from 1986 to 2015 (fig. 36). Divided into three 

distinct habitat states, 719 ha of predicted habitat were lost (55.0%), 25 ha were gained (1.9%), and 588 

ha (45.0%) persisted (fig. 37). Divided into nine channel zones, habitat losses varied from 18.79 percent

(zone 4) to 96.08 percent (zone 7) (fig. 38; table 6). Predicted changes in habitat by zone can be 

misleading without examining the amount of habitat involved. Specifically, zones 1–5 contained from

162.6 to 396.2 ha, compared with zones 6–9, which contained from 4.1 to 29.2 ha, making the 

magnitude of habitat loss much more significant in the first five zones.
62 

Figure 35. A habitat time series of predicted flycatcher habitat for upper Gila River, Arizona, from 1986 to 2015, 

plus the area of predicted flycatcher habitat after simulating a tamarisk leaf beetle infestation in 2015 (2015b). See 

table A5 for metadata on imagery used in this analysis. A satellite model was used to predict flycatcher habitat at a 

40-percent probability threshold.
63 

Figure 36. Map showing changes in predicted flycatcher habitat along the upper Gila River after conducting a 

tamarisk leaf beetle-impact simulation. The upper Gila River was divided into nine zones in order to quantify 

changes (see table 6 for the start and end points of each zone). A satellite model was used to predict flycatcher 

habitat at a 40-percent probability threshold.
64 

Figure 37. Graph showing three states of flycatcher habitat produced by a tamarisk leaf beetle-impact simulation 

along the upper Gila River (see fig. 10). Broken into three distinct habitat-change classes, 719 ha of habitat were 

lost, 25 ha were gained, while 588 ha persisted. A satellite model was used to predict flycatcher habitat at a 40-

percent probability threshold.

Figure 38. Graph showing simulated changes (percent) in predicted flycatcher habitat within 9 zones along the 

upper Gila River. A satellite model was used to predict flycatcher habitat at a 40-percent probability threshold.
65 

Table 6. Area of predicted flycatcher habitat before and after a tamarisk leaf beetle-impact simulation was done for 

the upper Gila River, Arizona.

[The satellite model produced the habitat estimates at a 40-percent probability threshold. The nine zones are delineated by 

rivers and creeks that are stored in the U.S. Geological Survey NHDPlus dataset. Zone: There are nine upper Gila River 

zones, as shown in figure 36. Habitat_pre: Flycatcher habitat (hectares) before leaf beetle-impact simulation was done.

Habitat_post: Flycatcher habitat (hectares) after leaf beetle-impact simulation was done. Change: Simulated change in 

predicted flycatcher habitat after a leaf beetle-impact simulation was done]

Zone

Habitat area (hectares) Change Description Habitat_pre Habitat_post

1 248.8 14.6 -94.14 Goodwin Wash to Black Rock Wash

2 396.2 203.1 -48.73 Black Rock Wash to Hot Springs Wash

3 264.5 110.9 -58.08 Hot Springs Wash to Markham Wash

4 187.7 152.5 -18.79 Markham Wash to Pack Wash

5 162.6 104.0 -36.03 Peck Wash to Butler Wash

6 29.2 20.7 -29.01 Butler Wash to Peterson Wash

7 4.6 0.2 -96.08 Peterson Wash to San Simon River

8 10.0 5.8 -42.34 San Simon River to Yuma Wash

9 4.1 1.5 -62.22 Yuma Wash to Bonita Creek

Discussion

Key Findings and a Path Forward

This report presents the first predictive, regionwide map of flycatcher breeding habitat, and there

are numerous findings likely to be of interest to natural resource managers, habitat modelers, and avian 

ecologists: 

1. A rangewide change-detection map that provides new insights into the formation and loss of 

predicted flycatcher breeding habitat; 

2. A hierarchical ranking of predicted flycatcher habitat at multiple scales—rangewide, 

statewide, FWS management units, critical-habitat reaches, and U.S. Geological Survey 7.5-

minute quadrangles; 

3. A set of spatially explicit maps that show predicted flycatcher habitat in binary or fiveprobability class formats for 2013–15; 

4. A rangewide accuracy assessment with a robust set of territories obtained from 2014; and

5. An assessment of beetle impacts along the lower Virgin River from 2010 to 2015, and

beetle-impact simulations along the lower Colorado and upper Gila Rivers in 2015. 
66 

Habitat Change Detection

Predicted flycatcher habitat increased 58.8 percent along a lateral west-east gradient from 2013 

to 2015, suggesting that regional weather patterns were responsible for greening and graying of riparian 

zones. The greening of riparian zones, as indicated by an increase in predicted flycatcher habitat from

2013 to 2015, was especially pronounced along the Rio Grande in southwestern Texas, the Pecos River 

in southeastern New Mexico, and the Rio Grande in western New Mexico (fig. 23). This pattern 

changed abruptly in Arizona where there was considerable variability in the state of flycatcher habitat, 

with some rivers exhibiting a checkerboard pattern of graying and greening from headwaters to mouth. 

The major tributaries of the Gila River (for example, San Pedro, Verde, San Francisco Rivers) all 

showed a considerable decrease in predicted flycatcher habitat (riparian areas became less lush), 

whereas the main stem Gila River showed increased lushness and flycatcher habitat along many 

reaches. There was another transition in California where most rivers and streams showed a pronounced 

decrease in predicted flycatcher habitat that indicated a loss of vigor in the riparian vegetation related to 

drought. These observations are strongly supported by U.S. Drought Monitor Maps (U.S. Department of 

Agriculture, 2015) that clearly show how the drought increased in the Western United States and 

decreased in the Eastern United States from 2013 to 2015 (fig. 39). The strong response of the satellite 

model to drought is not surprising because NDVI is a key component of the flycatcher model and is 

very sensitive to moisture and temperature (Williams and others, 2013). At finer scales, large interannual changes in predicted flycatcher habitat have been documented owing to drought, flooding, and 

surface and groundwater conditions (Hatten and others, 2010; Orr and others, 2014). 

The largest changes in predicted flycatcher habitat in 2015 were at reservoir inflows and riverine 

reaches. This is not surprising because fluctuating water levels caused by seasonal runoff or reservoir 

operations can result in optimal seed deposition and germination zones for riparian vegetation (Shafroth 

and others, 1998; Paradzick and Hatten, 2004). Young, lush vegetation in reservoir drawdown zones has

been shown to be heavily favored by breeding flycatchers (Paxton and others, 2007). At Roosevelt Lake 

in south-central Arizona, the satellite model provided habitat estimates that explained 79 percent of the 

inter-annual variability in a fluctuating flycatcher population from 1996 to 2004 (Hatten and others, 

2010). Additionally, flycatcher populations along Alamo Lake in eastern Arizona occur in areas 

predicted to be suitable by the satellite model (Hatten and Paradzick, 2003). In addition to lake deltas, 

riverine areas that have interval flooding can result in favorable flycatcher breeding habitat, such as the 

lower Bill Williams River or Rio Grande (Dockens and others, 2004; Hatten and Sogge, 2007).
67 

Figure 39. Drought maps portray a drought that intensifies to the west and decreases in the eastern portion of the 

study area from 2013 to 2015 (U.S. Department of Agriculture, 2015).
68 

Rangewide Ranking of Flycatcher Breeding Habitat

An intriguing part of this study was the identification of potentially new flycatcher breeding 

locations that have never been surveyed, or have not been visited in years. Given the scarcity of 

resources available to fund field surveys, this mapping effort could be used to identify and prioritize 

survey efforts. For example, the greening of the eastern part of flycatcher range could indicate that 

additional surveys are necessary along parts of the lower Rio Grande and Pecos River in eastern Texas. 

Although these areas only had only one documented breeding site (Durst and others, 2008), the satellite 

model indicates that the probability of favorable habitat has increased substantially from 2013 to 2015. 

Managers can be confident that the probability classes of the satellite model are a useful tool for 

prioritizing survey and restoration/enhancement efforts. For example, the likelihood of flycatcher 

territories is much greater in higher probability classes, with less than 12 percent in classes 1 or 2, thus, 

higher probability classes could be prioritized when survey funds are scarce (fig. 17). Using similar 

logic, it may be more economical to restore/enhance locations that are in the middle-probability range, 

such as class 3, because it is already much more likely to contain a flycatcher territory than classes 1 or 

2 (2.7 to 6 times, fig. 17). Such a strategy currently (2016) is being used along the upper Gila River in 

southeastern Arizona (near Fort Thomas), where tamarisk is being selectively removed and replaced 

with native willows before beetles arrive (Johnson and Calvo, 2014; Orr and others, 2014; Johnson, 

2016). 

Accuracy Assessment

A key question of this study was whether the performance of the satellite model was acceptable 

throughout the range of the flycatcher. The answer is “yes.” Specifically, the accuracy of the satellite 

model in 2014 was 88 percent (n=758) at a 40-percent cutpoint, with the remaining territories occurring 

on average only 74 m from predicted habitat. This accuracy assessment had a sample size that was 58

percent as large as the total known territories (1,299) rangewide since 1993 (Durst and others, 2008). 

Furthermore, these accuracy results compare favorably to other flycatcher mapping efforts in southcentral Arizona (Paxton and others, 2007; Hatten and others, 2010) and the Rio Grande, New Mexico 

(Hatten and Sogge, 2007). The fit of the model was respectable, with a clear exponential relationship 

between territory densities and five probability classes. 

A discussion about model accuracy is not complete without consideration of map accuracy 

because they are not the same thing (see Story and Congalton, 1986). For example, if the satellite model 

produces a 90-percent probability at a specific location (30 × 30-m cell), that does not mean that there is 

a 90-percent chance that a flycatcher territory will occur there, nor does it mean that it is flycatcher 

habitat. The only way to determine map accuracy from a user perspective is to randomly survey a set of 

locations in a given reach to determine their suitability (although there is subjectivity in visual 

classifications). After a given number of locations are randomly visited in each probability class, field 

results can be compared with model predictions (producer accuracy) using a classification table (error 

matrix). This approach was attempted after the first statewide mapping effort in 2001 (see Dockens and 

others, 2004), but it was limited because of budget constraints and logistical considerations (private 

access rights, etc.). However, results of Dockens and others (2004) are consistent with model accuracy 

because field-map misclassification decreased in high probability classes and increased in low 

probability classes. 
69 

The difficulties encountered in this field verification exercise actually reinforced the strengths of the 

satellite mode; mainly, that it can be consistently and cheaply applied without bias, provided the 

underlying ruleset (as described in the section, “Methods”) is not altered during model application. 

Beetle-Impact Analysis and Simulations

This is the first time that the flycatcher satellite model has been used to measure the impacts of 

the beetle on its habitat. The results are dramatic, largely because most of the lower Virgin River is 

comprised of tamarisk habitat. In areas that have reduced amounts of tamarisk and dominant native 

vegetation, losses to predicted flycatcher habitat likely will be less extreme. These results mirror 

previous studies that have shown large decreases in the vigor and lushness of tamarisk following beetle 

colonization in the upper Colorado and Dolores Rivers, Utah (Dennison and others, 2009), Humboldt 

and Walker Rivers, Nevada (Pattison and others, 2011), and lower Virgin River, Nevada and Arizona 

(Bateman and others, 2013; Nagler and others, 2014; Bateman and Johnson, 2015). Although these 

studies used different methods and satellite sensors, the conclusions of each study are very similar with

regard to the impacts of the beetle on tamarisk. This study went a step further by linking decreases in 

tamarisk vitality to flycatcher breeding habitat. 

The beetle-impact simulations provide information to managers of potential decreases in 

flycatcher habitat should beetles colonize the Lower Colorado or upper Gila Rivers. The best vegetation 

maps available were used for these reaches, but neither map had an accuracy assessment (Bangle and 

others, 2014; Neale and others, 2014), so caution is warranted when evaluating the beetle-impact

simulations. On the lower Virgin River, which was not a simulation, only 5.8 percent of predicted

flycatcher habitat persisted from 2010 to 2015 (beetles introduced in 2008), whereas simulation results 

indicated that 45 percent will persist along the upper Gila River and 64.1 percent will persist along the 

Lower Colorado River. On the upper Gila River, most of the predicted habitat was in zones 1–5 (fig. 

36), where predicted decreases in flycatcher habitat ranged from 18.8 to 94.1 percent. Beetles are not yet 

present on the upper Gila River; subtropical tamarisk beetles are currently (2016) predicted to arrive from 

western Texas and northern Arizona and New Mexico in 2–3 years (J. Tracy, Texas A&M University, 

written commun., April 2016) (fig. 2). The situation seems worse on the lower Colorado River, where 

northern tamarisk beetles have been moving steadily downstream and subtropical tamarisk beetles moving 

west will likely colonize the entire reach in the next 3–5 years (J. Tracy, Texas A&M University, written

commun., April 2016).

The beetle-impact simulations can help guide restoration/enhancement actions across the range 

of the flycatcher. Proactive efforts are occurring to minimize beetle effects on flycatcher habitat at 

multiple locations. The Gila Watershed Partnership is working proactively in zones 1–5 (fig. 36) by

initiating habitat restoration at sites that currently provide breeding habitat for Southwestern Willow 

Flycatchers and yellow-billed cuckoos (Coccyzus americanus). However, because of the anticipated 

arrival of the tamarisk leaf beetle (Diorhabda spp.), this organization is taking action now that will 

restore native habitat, minimize impacts of the beetle, and conserve this vulnerable riparian ecosystem

(Orr and others, 2014). A similar proactive effort is occurring 
70 

along the lower Colorado River, where Reclamation is replacing agricultural plantations with native 

vegetation more suitable to flycatchers and other native birds, such as yellow-billed cuckoos. The 

satellite model identified multiple native plantations as high-probability flycatcher habitat (fig. 32), 

providing support so that restoration and enhancement efforts can help minimize projected impacts from 

invasive beetles or natural disasters such as flooding, fire, or drought. 

Next Steps

Occupancy Modeling

Linking an occupancy model to the satellite model might significantly improve the ability to 

identify where flycatchers occur on an annual basis (MacKenzie and others, 2006). Satellite model 

simulations identify predicted flycatcher habitat and produce territory density estimates within five

probability classes, but they cannot estimate the likelihood that a flycatcher will occur. Large areas of 

predicted flycatcher habitat are vacant in some years but become occupied in subsequent years (U.S. 

Fish and Wildlife Service, 2002). Researchers have gathered considerable amounts of information 

related to flycatcher movement and habitat preferences (Paxton and others, 2007), but we lack the 

ability to predict which patches will become occupied as local or regional conditions change. The 

satellite model provides a foundation upon which an occupancy model could be built; the missing piece 

of the puzzle is a more complete accounting of flycatcher territories rangewide. I used 758 territories 

obtained in 2014 that were in GIS or spreadsheet format, but there were other locations handwritten on 

hard-copy maps that were not included in the analysis, indicating the importance of updating the 

regionwide flycatcher database on an annual basis (Durst and others, 2008). Three to five years of 

survey efforts coupled with satellite model results should provide sufficient data to produce an 

occupancy model (MacKenzie and others, 2006). A new suite of predictor variables would need to be 

created from the satellite model output, such as quantity of predicted habitat within different radii 

(Hatten and Paradzick, 2003), proximity to different sized creeks and rivers obtained from the NHDPlus 

database (McKay and others, 2012), location (FWS management unit, 7.5-minute quadrangle, HUCHydrologic Unit Code, easting), and number of years a site has been above or below a given probability 

threshold (Hatten and others, 2010). 

Automated Cloud-Based Mapping

The methods and results presented in this report could be used in the development of an 

automated cloud-based mapping application. Specifically, the satellite model algorithms could be 

applied automatically to Landsat scenes that are being collected monthly by Google EarthTM mapping 

service from the U.S. Geological Survey EROS Data Center (Google EarthTM, 2016a). Created in 

cooperation with the U.S. Geological Survey, Google EarthTM Engine was developed as a planetary 

application that could track changes in forest cover (Regalado, 2010). Through Google Earth Engine, 

users could access flycatcher habitat maps with an interactive Internet mapping application, similar to 

the recent Google EarthTM Engine application developed for sage grouse (U.S. Department of 

Agriculture, 2016), specifying an area of interest (for example, polygon) and date (could be multiple 

years). Google EarthTM Engine could enable quick and efficient rangewide maps to be created and the 

development of habitat time series for any locations from 1986 to the present (Google EarthTM, 2016b). 
71 

Acknowledgments

I thank Mark Sogge (U.S. Geological Survey) and Carol Evans (U.S. Bureau of Reclamation) 

for helping with project initiation and coordination, and Thomas Batt (USGS) for downloading and 

organizing Landsat data. I thank Susan Sferra, Greg Beatty, and Scott Durst of the U.S. Fish and 

Wildlife Service for providing logistical, technical, and administrative support. I also thank the 

following individuals for providing verification data: Sonja Kokos (Bureau of Reclamation), Vicki Ryan 

(U.S. Fish and Wildlife Service), Darrell Ahlers (Bureau of Reclamation), Barbara Kus and Scarlett 

Howell (USGS), Matthew Johnson (Northern Arizona University), Lacey Green (California Department 

of Fish and Wildlife), and Martha Schuman (The Nature Conservancy). I also thank Eben Paxton, Linda 

Rogers, John Osias, and Barbara Ralston (U.S. Geological Survey), and Matthew Johnson (Northern 

Arizona University), for providing useful comments that greatly improved this paper, and Luis Menoyo 

(U.S. Geological Survey) for graphical support. A hearty thanks to all the nameless field technicians that 

labored to collect territory and nest locations throughout the flycatcher’s range. Lastly, I thank James

Tracy (Texas A&M University) for providing useful information related to tamarisk beetles. 

References Cited

Anderson, B.W., and Ohmart, R.D., 1984, Lower Colorado River riparian methods of quantifying vegetation 

13 communities to prepare type maps: Final report prepared for the Bureau of Reclamation, Lower 

Colorado Region, Boulder City, Nevada.

Bateman, H.L., Nagler, P.L., and Glenn, E.P., 2013, Plot- and landscape-level changes in climate and 

vegetation following defoliation of exotic saltcedar (Tamarix sp.) from the biocontrol agent 

Diorhabda carinulata along a stream in the Mojave Desert (USA): Journal of Arid Environments, v. 

89, p. 16–20. 

Bateman, H.L., and Johnson, M.J., 2015, Effects of biocontrol and restoration on wildlife in 

Southwestern riparian habitats: Report submitted to Southern Rockies and Desert Landscape 

Conservation Cooperative. 

Bangle, D., Ahlbrandt, R., Baker, M., and Milliken, J., 2014, 2010–2011 vegetation mapping along the 

lower Colorado River . Bureau of Reclamation, unpub. report, 74 p.

Bean, D.W., Dudley, T., and Huletine, K., 2013, Bring on the beetles!—The history and impact of 

tamarisk biological control, in Sher, A., and Quigley, M.F., eds., Tamarix—A case study of ecological 

change in the American West: New York, Oxford University Press, p. 377–403. 

Chander, G., Markham, B.L., and Helder, D.L., 2009, Summary of current radiometric calibration 

coefficients for Landsat MSS, TM, ETM+, and EO-1 sensors: Remote Sensing of Environment, v. 

113, p. 893–903. 

Chew, M.K., 2013, Tamarisk introduction, naturalization, and control in the United States, 1818–1952,

in Sher, A., and Quigley, M.F., eds., Tamarix—A case study of ecological change in the American 

West; New York, Oxford University Press, p. 269–286. 

DeLoach, C.J., Carruthers, R., 2004, Biological control programs for integrated invasive plant 

management, in Proceedings of Weed Science Society of America Meeting, Kansas City, Missouri: 

Weed Science Society of America (CD-ROM), 17 p. 

Dennison, P.E, Nagler, P.L., Hultine, K.R., Glenn, E.P., and Ehleringer, J.R., 2009, Remote monitoring 

of tamarisk defoliation and evapotranspiration following saltcedar leaf beetle attack: Remote Sensing 

of Environment, v. 113, p. 1,462–1,472. 

Dobbs, R. C., Huizinga, M., Edwards, C.N., and Fridell, R., 2012, Status, reproductive success, and 

habitat use of Southwestern willow flycatchers on the Virgin River, Utah, 2008–2011: Publication 

Number 12-36, Utah Division of Wildlife Resources. 
72 

Dockens, P.E.T., Paradzick, C.E., and Hatten, J.R., 2004, Application of a Southwestern willow 

flycatcher GIS-based habitat model—An estimate of breeding habitat in Arizona, 2001, in Dockens, 

P.E.T., and Paradzick, C.E., eds., Mapping and monitoring Southwestern willow flycatcher breeding 

habitat in Arizona—A remote sensing approach: Phoenix, Arizona Game and Fish Department 

Nongame and Endangered Wildlife Program Technical Report 223, p. 28–59. 

Durst, S.L., Sogge, M.K., Stump, S.D., Walker, H,A, Kus, B.E., and Sferra, S.J., 2008. Southwestern 

willow flycatcher breeding sites and territory summary—2007: U.S. Geological Survey Open-File 

Report 2008-1303, http://pubs.usgs.gov/of/2008/1303. 

Google Earth, 2016a, A planetary-scale platform for Earth science data and analysis: Google Earth 

Engine Web site, accessed February 29, 2016, at https://earthengine.google.com/. 

Google Earth, 2016b, Timelapse: Google Earth Engine Web site, accessed February 29, 2016, at 

https://earthengine.google.com/timelapse/. 

Hatten, J.R., 2014, Mapping and monitoring Mount Graham red squirrel habitat with lidar and Landsat 

imagery: Ecological Modelling, v. 289, p. 106–123. 

Hatten, J.R., and Paradzick, C.E., 2003, A multiscaled model of southwestern willow flycatcher 

breeding habitat: Journal of Wildlife Management, v. 67, p. 774–788. 

Hatten, J.R., Paxton, E.H., and Sogge, M.K., 2010, Modeling the dynamic habitat and breeding 

population of Southwestern willow flycatcher: Ecological Modelling, v. 221, p. 1,674–1,686. 

Hatten, J.R., and Sogge, M.K., 2007, Using a remote sensing/GIS model to predict Southwestern willow 

flycatcher breeding habitat along the Rio Grande, New Mexico: U.S. Geological Survey Open-File 

Report 2007-1207, http://pubs.usgs.gov/of/2007/1207/. 

Hultine, K.R., Belnap, J., van Riper, C., III, Ehleringer, J.R., Dennison, P.E., Lee, M.E., Nagler, P.E., 

Snyder, K.A., Uselman, S.M., and West, J.B., 2009, Tamarisk biocontrol in the Western United 

States—Ecological and societal implications: Frontiers in Ecology and the Environment e-View, doi: 

10.1890/090031. 

Hultine, K.R., Dudley, T.L., Koepke, D.F., Bean, D.W., Glenn, E.P., and Lambert, A.M., 2015, Patterns 

of herbivory-induced mortality of a dominant non-native tree/shrub (Tamarix spp.) in a Southwestern 

US watershed: Biological Invasions, doi: 10.1007/s10530-014-0829-4. 

Jamison, L.R., van Riper, C., III, and Bean, D.W., 2015, The influence of Tamarix ramosissima

defoliation on population movements of the Northern tamarisk beetle (Diorhabda carinulata) within 

the Colorado Plateau, in Huenneke, L.F., van Riper, C., III, and Hays-Gilpin, K.A., The Colorado 

Plateau VI—Science and management at the landscape scale: Tucson, University of Arizona Press, p. 

289–300, Project MUSE, https://muse.jhu.edu/, 

Johnson, M.J., 2016, Gila River, AZ Southwestern willow flycatcher and yellow-billed cuckoo surveys 

and nest monitoring for proposed restoration sites, 2015: Annual report submitted to Gila Watershed 

Partnership of Arizona and U.S. Fish and Wildlife Service,. 42 p. 

Johnson, M.J., and Calvo, C.M., 2014, Virgin River tamarisk removal and restoration and its effects on 

Southwestern willow flycatchers—Virgin River, Arizona assessment and recommendations: Annual 

report submitted to Arizona State Forestry Division, 113 p. 

Johnson, M.J., and Jamison, L.R.. 2015, Tamarisk leaf beetle species distribution in Rio Grande River 

watershed, New Mexico, 2015: Bosque Ecosystem Monitoring Program, Annual Report submitted to 

the Army Corp of Engineers, Albuquerque, New Mexico, 71 p. 

Lowry, J.H, Jr., Ramsey, R.D., Boykin, K., Bradford, D., Comer, P., Falzarano, S., Kepner, W., Kirby, 

J., Langs, L., Prior-Magee, J., Manis, G., O’Brien, L., Sajwaj, T., Thomas, K.A., Rieth, W., Schrader, 

S., Schrupp, D., Schulz, K., Thompson, B., Velasquez, C., Wallace, C., Waller, E., and Wolk, B., 

2005, Southwest regional gap analysis project—Final report on land cover mapping methods: Logan, 

Utah State University, RS/GIS Laboratory. 
73 

MacKenzie, D.I., Nichols, J.D., Royle, J.A., Pollock, K.H., Bailey, L.L., Hines, J.E., 2006, Occupancy 

estimation and modeling—Inferring patterns and dynamics of species occurrence: San Diego, 

Academic Press.

McKay, L., Bondelid, T., Rea, A., Johnston, C., Moore, R., and Dewald, T., 2012, NHDPlus Version 

2—User Guide: Prepared for U.S. Environmental Protection Agency, 159 p. 

Nagler, P.L., Pearlstein, S., Glenn, E.P., Brown, T.B., Bateman, H.L., Bean, D.W., and Hultine, K.R., 

2014, Rapid dispersal of saltcedar (Tamarix spp.) biocontrol beetles (Diorhabda carinulata) on a 

desert river detected by phenocams, MODIS imagery and ground observations: Remote Sensing of 

Environment, v. 140, p. 206–219. 

Neale, C.M.U., Akasheh, O., and McGinty, C., 2014, Mapping the geomorphology and vegetation in the 

Gila River riparian corridor to support restoration planning and assessment: Logan, Utah State 

University, Department of Civil and Environmental Engineering,, unpub. report, 17 p. 

Orr, B.K., Leverich, G.T., Diggory, Z.E., Dudley, T.L., Hatten, J.R., Hultine, K.R., Johnson, M.P., and 

Orr, D.A., 2014, Riparian restoration framework for the upper Gila River, Arizona: Prepared for the 

Gila Watershed Partnership of Arizona, 57 p. with appendixes. 

Paradzick, C.E., and Hatten, J.R., 2004, Southwestern willow flycatcher habitat change detection

analysis—Roosevelt Lake, Arizona, 1985–2001 in Dockens, P.E.T., and Paradzick, C.E., eds., 

Mapping and monitoring Southwestern willow flycatcher habitat in Arizona—A remote sensing 

approach: Phoenix, Arizona Game and Fish Department Nongame and Endangered Wildlife Program 

Technical Report 223, p. 60–82. 

Pattison, R.R., D’Antonio, C.M., Dudley, T.L., Allander, K.K., and Rice, R., 2011, Early impacts of 

biological control on canopy cover and water use of the invasive saltcedar tree (Tamarix spp.) in 

western Nevada, USA: Oecologia, v. 165, p. 605–616. 

Paxton, E.H., Sogge, M.K., Durst, S.L., Theimer, T.C., and Hatten, J.R., 2007, Spatial Modeling, chap. 

6 of The ecology of the Southwestern willow flycatcher in central Arizona—A 10-year synthesis 

report: U.S. Geological Survey Open-File Report 2007-1381, p 97–124. 

Paxton, E.H., Theimer, T.C., and Sogge , M.K., 2011, Tamarisk biocontrol using tamarisk beetles—

Potential consequences for riparian birds in the Southwestern United States: The Condor, v. 113, no, 

2, p. 255–265. 

Regalado, A., 2010, New Google Earth Engine: Science magazine online, accessed February 29, 2016, 

at http://www.sciencemag.org/news/2010/12/new-google-earth-engine. 

Shafroth, P.B., Auble, G.T., Stromberg, J.C., and Patten, D.T., 1998, Establishment of woody riparian 

vegetation in relation to annual patterns of stream flow, Bill Williams River, Arizona: Wetlands, v. 

18, p. 577–590. 

Sogge, M.K., Ahlers, D., and Sferra, S.J., 2010, A natural history summary and survey protocol for the 

Southwestern willow flycatcher: U.S. Geological Survey Techniques and Methods, book 2, chap. 

A10, http://pubs.usgs.gov/tm/tm2a10/. 
74 

Sogge, M.K., and Marshall, R.M., 2000, Chapter 5: A survey of current breeding habitats, in Finch, 

D.M., and Stoleson, S.H., eds., Status, ecology, and conservation of the Southwestern Willow 

Flycatcher: U.S. Forest Service Rocky Mountain Research Station General Technical Report-60, p. 

43-56.

Sogge, M.K., Paxton, E.H., and van Riper, C., III, 2013, Tamarisk in riparian woodlands—A bird’s eye 

view, in Sher A., and Quigley, M.F., eds., Tamarix—A case study of ecological change in the 

American West: New York, Oxford University Press, p. 189–206. 

Sogge, M.K., Sferra, S.J., and Paxton, E.H., 2008, Tamarix as habitat for birds—Implications for 

riparian restoration in the Southwestern United States: Restoration Ecology, v. 16, p. 146–154. 

Story, M., Congalton, R.G., 1986. Accuracy assessment: a user’s perspective. Photogramm. Eng. 

Remote Sens. 52, 397–399. 

U.S. Fish and Wildlife Service, 1995, Final rule determining endangered status for the Southwestern 

willow flycatcher: Federal Register, v. 60, no. 38, p. 10694–10715. 

U.S. Fish and Wildlife Service, 2002, Southwestern willow flycatcher final recovery plan: U.S. Fish and 

Wildlife Service, Albuquerque, New Mexico. 

U.S. Fish and Wildlife Service, 2013, Designation of critical habitat for Southwestern willow 

flycatcher—Final rule: Federal Register, v. 78, no. 2., p. 344–534. [Also available at

http://www.regulations.gov, Docket No. FWS–R2–ES–2011–0053.] 

U.S. Department of Agriculture, 2015, U.S. Drought Monitor (home page): Lincoln, Nebraska, National 

Drought Mitigation Center Web site, accessed February 24, 2016, at http://droughtmonitor.unl.edu/.

U.S. Department of Agriculture, 2016, Sage grouse initiative—Interactive viewer: Sage Grouse 

Initiative Web site, accessed February 29, 2016, at http://map.sagegrouseinitiative.com/. 

U.S. Geological Survey, 2005, Global Visualization (GloVis) viewer: U.S. Geological Survey Fact 

Sheet 2005-3118, 1 p., accessed February 11, 2016, at http://glovis.usgs.gov/. 

Williams, A.P., Allen, D.A., Macalady, A.K., Griffin, D., Woodhouse, C.A., Meko, D.M., and others, 

2013, Temperature as a potent driver of regional forest drought stress and tree mortality: Nature 

Climate Change, v. 3, p. 292–297.
75 

Appendix A. Metadata for Landsat Scenes Used in Regionwide Mapping 

and Habitat Time Series

Table A1. Metadata for 57 Landsat scenes used to model flycatcher breeding habitat rangewide, 

Southwestern United States, 2013.

[UTM zone: Universal Transverse Mercator zone]

Date Scene_ID UTM zone Clouds Quality

ranking

6/1/2013 LC80300372013152LGN00 14 0.12 9

6/1/2013 LC80300382013152LGN00 13 3.88 9

6/1/2013 LC80300392013152LGN00 13 0.53 9

6/1/2013 LC80300402013152LGN00 13 0.23 9

6/8/2013 LC80310362013159LGN00 13 0.06 9

6/8/2013 LC80310372013159LGN00 13 0.05 9

6/8/2013 LC80310382013159LGN00 13 0.12 9

46/8/2013 LC80310392013159LGN00 13 0.05 9

6/8/2013 LC80310402013159LGN00 13 0.01 9

5/14/2013 LC80320352013134LGN03 13 0.24 9

6/15/2013 LC80320362013166LGN00 13 0.2 9

6/15/2013 LC80320372013166LGN00 13 1.7 9

6/15/2013 LC80320382013166LGN00 13 0.3 9

6/15/2013 LC80320392013166LGN00 13 0.35 9

6/6/2013 LC80330342013157LGN00 13 2.15 9

6/22/2013 LC80330352013173LGN00 13 3.02 9

6/22/2013 LC80330362013173LGN00 13 2.06 9

6/6/2013 LC80330372013157LGN00 13 5.15 9

6/6/2013 LC80330382013157LGN00 13 0.09 9

7/31/2013 LC80340342013212LGN00 13 9.25 9

6/29/2013 LC80340352013180LGN00 13 5.11 9

6/29/2013 LC80340362013180LGN00 13 7.35 9

5/28/2013 LC80340372013148LGN00 13 2.1 9

6/13/2013 LC80340382013164LGN00 12 0.07 9

6/20/2013 LC80350342013171LGN00 12 0.06 9

6/20/2013 LC80350352013171LGN00 12 0.01 9

6/20/2013 LC80350362013171LGN00 12 0.02 9

6/20/2013 LC80350372013171LGN00 12 0.01 9

6/20/2013 LC80350382013171LGN00 12 0.04 9

6/11/2013 LC80360342013162LGN00 12 0 9

6/27/2013 LC80360352013178LGN01 12 0.01 9

6/27/2013 LC80360362013178LGN01 12 0.02 9

6/27/2013 LC80360372013178LGN01 12 0.11 9
76 

Date Scene_ID UTM zone Clouds Quality

ranking

6/11/2013 LC80360382013162LGN00 12 0.18 9

6/18/2013 LC80370342013169LGN01 12 0.02 9

6/18/2013 LC80370352013169LGN01 12 0 9

6/18/2013 LC80370362013169LGN01 12 0.02 9

6/18/2013 LC80370372013169LGN01 12 0.24 9

6/18/2013 LC80370382013169LGN01 12 0.13 9

6/25/2013 LC80380342013176LGN00 12 0.15 9

6/9/2013 LC80380352013160LGN00 12 0.01 9

6/9/2013 LC80380362013160LGN00 12 0.06 9

6/25/2013 LC80380372013176LGN00 11 0.05 9

6/25/2013 LC80380382013176LGN00 11 5.93 9

6/16/2013 LC80390342013167LGN00 11 0.02 9

6/16/2013 LC80390352013167LGN00 11 0.26 9

6/16/2013 LC80390362013167LGN00 11 0.08 9

6/16/2013 LC80390372013167LGN00 11 0.34 9

6/16/2013 LC80390382013167LGN00 11 6.23 9

6/23/2013 LC80400342013174LGN00 11 0.12 9

6/23/2013 LC80400352013174LGN00 11 0.63 9

6/23/2013 LC80400362013174LGN00 11 0.96 9

5/22/2013 LC80400372013142LGN01 11 27.79 9

6/14/2013 LC80410342013165LGN00 11 1 9

6/14/2013 LC80410352013165LGN00 11 1.05 9

5/29/2013 LC80410362013149LGN00 11 2.77 9

5/29/2013 LC80410372013149LGN00 11 7.66 9

6/21/2013 LC80420342013172LGN00 11 1.23 9

6/21/2013 LC80420352013172LGN00 11 0.18 9

6/21/2013 LC80420362013172LGN00 11 1.5 9

6/12/2013 LC80430342013163LGN00 10 0.49 9

6/28/2013 LC80430352013179LGN00 10 11.55 9

6/19/2013 LC80440342013170LGN00 10 0.47 9
77 

Table A2. Metadata for 57 Landsat scenes used to model flycatcher breeding habitat rangewide, 

Southwestern United States, 2014. 

[UTM zone: Universal Transverse Mercator zone] 

Date Scene_ID UTM

zone

Clouds Quality

ranking

6/4/2014 LC80300372014155LGN00 14 0.04 9

6/4/2014 LC80300382014155LGN00 13 0.04 9

6/4/2014 LC80300392014155LGN00 13 0.02 9

6/4/2014 LC80300402014155LGN00 13 0.02 9

6/11/2014 LC80310362014162LGN00 13 4.11 9

5/10/2014 LC80310372014130LGN00 13 0.05 9

6/11/2014 LC80310382014162LGN00 13 4.39 9

6/11/2014 LC80310392014162LGN00 13 0.09 9

6/27/2014 LC80310402014178LGN00 13 0.26 9

5/17/2014 LC80320342014137LGN00 13 0.20 9

5/17/2014 LC80320352014137LGN00 13 0.01 9

6/2/2014 LC80320362014153LGN00 13 0.02 9

6/2/2014 LC80320372014153LGN00 13 0.46 9

6/2/2014 LC80320382014153LGN00 13 0.35 9

6/2/2014 LC80320392014153LGN00 13 0.00 9

6/9/2014 LC80330342014160LGN00 13 0.35 9

6/9/2014 LC80330352014160LGN00 13 0.04 9

6/9/2014 LC80330362014160LGN00 13 0.12 9

6/9/2014 LC80330372014160LGN00 13 2.49 9

6/9/2014 LC80330382014160LGN00 13 0.11 9

6/16/2014 LC80340342014167LGN00 13 0.87 9

6/16/2014 LC80340352014167LGN00 13 0.03 9

6/16/2014 LC80340362014167LGN00 13 3.60 9

5/31/2014 LC80340372014151LGN00 13 4.02 9

6/16/2014 LC80340382014167LGN00 12 4.07 9

6/23/2014 LC80350342014174LGN00 12 18.09 9

6/23/2014 LC80350352014174LGN00 12 0.61 9

6/7/2014 LC80350362014158LGN00 12 0.02 9

6/7/2014 LC80350372014158LGN00 12 0.01 9

6/23/2014 LC80350382014174LGN00 12 0.67 9

5/13/2014 LC80360342014133LGN00 12 1.62 9

6/14/2014 LC80360352014165LGN00 12 0.01 9

6/14/2014 LC80360362014165LGN00 12 0.03 9

6/14/2014 LC80360372014165LGN00 12 0.45 9

6/14/2014 LC80360382014165LGN00 12 0.34 9

6/5/2014 LC80370342014156LGN01 12 0.73 9

6/5/2014 LC80370352014156LGN01 12 0.01 9
78 

Date Scene_ID UTM

zone

Clouds Quality

ranking

6/21/2014 LC80370362014172LGN00 12 0.29 9

6/21/2014 LC80370372014172LGN00 12 0.28 9

6/21/2014 LC80370382014172LGN00 12 0.08 9

6/28/2014 LC80380342014179LGN00 12 0.06 9

6/28/2014 LC80380352014179LGN00 12 0.01 9

6/28/2014 LC80380362014179LGN00 12 0.06 9

6/28/2014 LC80380372014179LGN00 11 0.04 9

6/28/2014 LC80380382014179LGN00 11 1.79 9

6/19/2014 LC80390342014170LGN00 11 0.04 9

6/19/2014 LC80390352014170LGN00 11 0.26 9

6/19/2014 LC80390362014170LGN00 11 0.06 9

6/19/2014 LC80390372014170LGN00 11 0.32 9

6/19/2014 LC80390382014170LGN00 11 1.07 9

5/25/2014 LC80400342014145LGN00 11 0.36 9

6/10/2014 LC80400352014161LGN00 11 2.04 9

6/10/2014 LC80400362014161LGN00 11 0.41 9

5/25/2014 LC80400372014145LGN00 11 20.97 9

6/17/2014 LC80410342014168LGN00 11 2.11 9

6/17/2014 LC80410352014168LGN00 11 1.13 9

6/17/2014 LC80410362014168LGN00 11 0.93 9

6/17/2014 LC80410372014168LGN00 11 13.09 9

6/24/2014 LC80420342014175LGN00 11 1.25 9

6/24/2014 LC80420352014175LGN00 11 0.23 9

6/24/2014 LC80420362014175LGN00 11 15.83 9

6/15/2014 LC80430342014166LGN00 10 3.53 9

5/30/2014 LC80430352014150LGN00 10 9.00 9

6/22/2014 LC80440342014173LGN00 10 30.02 9
79 

Table A3. Metadata for 57 Landsat scenes used to model flycatcher breeding habitat rangewide, 

Southwestern United States, 2015.

[UTM zone: Universal Transverse Mercator zone]

Path/

row Scene_ID Acquisitiondate

UTM

zone

Clouds Quality

ranking

3038 LC80300382015206LGN00 7/25/2015 13 0.09 9

3039 LC80300392015206LGN00 7/25/2015 13 0.80 9

3040 LC80300402015158LGN00 6/07/2015 13 1.43 9

3136 LC80310362015181LGN00 6/30/2015 13 0.06 9

3137 LC80310372015165LGN00 6/14/2015 13 0.08 9

3138 LC80310382015165LGN00 6/14/2015 13 0.14 9

3139 LC80310392015197LGN00 7/16/2015 13 0.17 9

3140 LC80310402015197LGN00 7/16/2015 13 0.07 9

3235 LC80320352015172LGN00 6/21/2015 13 0.00 9

3236 LC80320362015172LGN00 6/21/2015 13 0.00 9

3237 LC80320372015172LGN00 6/21/2015 13 0.51 9

3238 LC80320382015172LGN00 6/21/2015 13 0.38 9

3239 LC80320392015172LGN00 6/21/2015 13 0.21 9

3334 LC80330342015227LGN00 8/15/2015 13 1.20 9

3335 LC80330352015227LGN00 8/15/2015 13 0.88 9

3336 LC80330362015227LGN00 8/15/2015 13 0.31 9

3337 LC80330372015227LGN00 8/15/2015 13 5.51 9

3338 LC80330382015163LGN00 6/12/2015 13 0.28 9

3434 LC80340342015154LGN00 6/03/2015 13 1.82 9

3435 LC80340352015170LGN00 6/19/2015 13 0.26 9

3436 LC80340362015170LGN00 6/19/2015 13 0.13 9

3437 LC80340372015170LGN00 6/19/2015 13 0.31 9

3438 LC80340382015170LGN00 6/19/2015 12 0.03 9

3534 LC80350342015177LGN00 6/26/2015 12 1.20 9

3535 LC80350352015177LGN00 6/26/2015 12 4.99 9

3536 LC80350362015177LGN00 6/26/2015 12 10.50 9

3537 LC80350372015177LGN00 6/26/2015 12 6.47 9

3538 LC80350382015177LGN00 6/26/2015 12 0.42 9

3634 LC80360342015168LGN00 6/17/2015 12 0.03 9

3635 LC80360352015168LGN00 6/17/2015 12 0.00 9

3636 LC80360362015168LGN00 6/17/2015 12 0.02 9

3637 LC80360372015168LGN00 6/17/2015 12 0.28 9

3638 LC80360382015152LGN00 6/1/2015 12 0.20 9

3734 LC80370342015159LGN00 6/8/2015 12 0.80 9

3735 LC80370352015191LGN00 7/10/2015 12 0.73 9

3736 LC80370362015175LGN00 6/24/2015 12 6.35 9

3737 LC80370372015191LGN00 7/10/2015 12 3.07 9

3738 LC80370382015191LGN00 7/10/2015 12 2.72 9
80 

Path/

row Scene_ID Acquisitiondate

UTM

zone

Clouds Quality

ranking

3834 LC80380342015198LGN00 7/17/2015 12 0.14 9

3835 LC80380352015182LGN00 70/1/2015 12 0.01 9

3836 LC80380362015198LGN00 7/17/2015 12 0.76 9

3837 LC80380372015166LGN00 6/15/2015 11 0.05 9

3838 LC80380382015166LGN00 6/15/2015 11 1.37 9

3934 LC80390342015173LGN00 6/22/2015 11 0.03 9

3935 LC80390352015173LGN00 6/22/2015 11 0.26 9

3936 LC80390362015173LGN00 6/22/2015 11 0.07 9

3937 LC80390372015173LGN00 6/22/2015 11 0.73 9

4034 LC80400342015196LGN00 7/15/2015 11 0.09 9

4035 LC80400352015196LGN00 7/15/2015 11 0.32 9

4036 LC80400362015196LGN00 7/15/2015 11 0.34 9

4037 LC80400372015196LGN00 7/15/2015 11 0.82 9

4134 LC80410342015171LGN00 6/20/2015 11 2.41 9

4135 LC80410352015171LGN00 6/20/2015 11 0.86 9

4136 LC80410362015171LGN00 6/20/2015 11 0.77 9

4137 LC80410372015171LGN00 6/20/2015 11 3.04 9

4234 LC80420342015194LGN00 7/13/2015 11 0.86 9

4235 LC80420352015194LGN00 7/13/2015 11 0.16 9

4236 LC80420362015226LGN00 8/14/2015 11 0.25 9
81 

Table A4. Metadata for habitat-time-series analysis of Virgin River, Nevada and Arizona.

[Satellite model was used to predict flycatcher habitat at a 40-percent probability threshold. For this analysis, two 

Landsat scenes were mosaicked together each year (paths/rows [PRs] 3934 and 3935)]

Year PR 3934

(Julian date)

Acquisition

date

PR 3935

(Julian date)

Acquisition

date

1986 221 8.8.86 205 7.23.86

1987 160 6.9.87 160 6.9.87

1988 163 6.11.88 195 7.13.88

1989 165 6.14.89 181 6.30.89

1990 184 7.2.90 184 7.2.90

1991 171 6.20.91 171 6.20.91

1992 182 6.30.92 182 6.30.92

1993 192 7.11.93 192 7.11.93

1994 179 6.27.94 179 6.27.94

1995 182 7.1.95 182 7.1.95

1996 169 6.17.96 169 6.17.96

1997 187 7.6.97 171 6.20.97

1998 174 6.23.98 158 7.7.98

1999 177 6.26.99 177 6.26.99

2000 180 6.28.00 148 5.27.00

2001 166 6.15.01 166 6.15.01

2002 169 6.18.02 169 6.18.02

2003 172 6.21.03 172 6.21.03

2004 159 6.7.04 159 6.7.04

2005 177 6.26.05 177 6.26.05

2006 180 6.29.06 180 6.29.06

2007 167 6.16.07 167 6.16.07

2008 170 6.18.08 170 6.18.08

2009 172 6.21.09 172 6.21.09

2010 175 6.24.10 175 6.24.10

2011 178 6.27.11 178 6.27.11

2013 167 6.16.13 167 6.16.13

2014 170 6.19.14 170 6.19.14

2015 173 6.22.15 173 6.22.15
82 

Table A5. Metadata for Landsat scenes used in habitat-time-series analysis, by year, upper Gila River, 

Arizona and California.

[All scenes are path 35 and row 37, UTM Zone 12. NA, not available]

Year Julian

day

Month and

day

1986 161 June 10

1987 180 June 28

1988 167 June 15

1989 169 June 18

1990 236 August 24

1991 159 June 8

1992 178 June 26

1993 180 June 29

1994 167 June 16

1995 170 June 19

1996 157 June 5

1997 175 June 24

1998 178 June 27

1999 181 June 30

2000 168 June 16

2001 170 June 19

2002 173 June 22

2003 176 June 25

2004 163 June 11

2005 165 June 14

2006 168 June 17

2007 123 May 3

2008 158 June 6

2009 128 May 8

2010 147 May 26

2011 166 June 15

2012 NA NA

2013 171 June 20

2014 158 June 7

2015 177 June 26
83 

Appendix B. Normalized Difference Vegetation Index (NDVI) Conversions 

for Landsat 8

Landsat 5

The flycatcher satellite model was developed with imagery acquired from the Landsat 5 

Thematic Mapper sensor that was in continuous operation from 1984 through 2011. No 

adjustments to Landsat 5 images are required before their use in the satellite model; just follow a 

few simple steps:

1. Preview and download the desired imagery from the EROS Data Center with the 

GLOVIS tool (http://glovis.usgs.gov/), 

2. calculate NDVI with bands 3 and 4, 

3. Create the three vegetation variables, and

4. Apply equation 1 (see section, “Methods”).

It is not advisable at this time to adjust Landsat 5 imagery to top-of-atmosphere (TOA) 

reflectance (Chander and others, 2009) because equation 1 will produce unreliable results. 

If one selects an image that was acquired with Landsat 8 (from late 2012 to present), the 

image must first be adjusted to TOA reflectance before calculating NDVI. Following TOA 

adjustment and NDVI calculation, a simple two-step regression is applied to the NDVI values to 

convert them to Landsat 5 NDVI values. Without applying the TOA adjustment or regressions, 

the satellite model will not produce results similar to those from 1984 to 2011. Some of the 

material presented in the following section (eq. B1) was published in a previous study (Hatten, 

2014). However, a significant alteration has been made (eq. B2 is new) in order for the flycatcher 

satellite model to run correctly on Landsat 8. Specifically, Landsat 5 uses bands 3 and 4 to 

calculate NDVI, whereas Landsat 8 uses bands 4 and 5. The bandwidth of Band 3 (red) for 

Landsat 5 is 0.63–0.69 μm, whereas band 4 (reflected IR) is 0.76–0.90 μm. In contrast, the 

Landsat 8 red band (band 4) is 0.64–0.67 μm, whereas its reflected IR band (band 5) is 0.85–0.88 

μm. These differences in bandwidth present problems when calculating NDVI because they 

produce slightly different NDVI values that affect model probabilities used to produce flycatcher 

habitat suitability maps. 

To successfully run the flycatcher satellite model after 2011 (Landsat 5 failed), the 

Landsat 8 NDVI values were calibrated against those of Landsat 5. Landsat 5 stopped working 

during autumn 2011, whereas Landsat 8 did not become operational until May 2013, making 

calibrating the two Landsat NDVI values impossible. A workaround was built by creating a 

statistical bridge from Landsat 8 to 7 NDVI values for 2013, and then from Landsat 7 to 5 NDVI

values for 2011. This was accomplished by first adjusting Landsat 8 and 7 to TOA reflectance 

(Chander and others, 2009), but not adjusting Landsat 5 to TOA (rationale provided below), and 

applying a two-step linear regression; step 1 converted Landsat 8 to Landsat 7 NDVI values, and

step two converted Landsat 7 to Landsat 5 NDVI values. Following these corrections, the three 

vegetation variables were created with the modified NDVI values and the flycatcher model was 

created by applying equation 1. Testing of these methods at multiple locations and times showed 

that they produced consistent results, whereas unadjusted imagery produced inconsistent results. 
84 

The NDVI adjustment procedure was not straightforward because Landsat 7 has an 

ongoing problem with one of its sensor arrays that produces no-data stripes in certain areas (this 

is why the use of Landsat 7 data is not recommended unless absolutely necessary). Every 

Landsat 7 scene has an error-free path (a large strip with no striping problems), so a stripe-free 

area was selected along the upper Gila River (fig. B1; Area 1) and a regression was developed 

between Landsat 8 and 7 NDVI values. I switched areas when I converted Landsat 7 to 5 NDVI 

values (fig. B1; Area 2) because of a better match on acquisition dates. 

After converting Landsat 7 and 8 images to TOA reflectance, NDVI was calculated for 

both images. Next, 5,000 random points were generated around the Salt River where there was 

much spectral variability because of desert and riparian areas. A geographic information system 

was used to attribute the 5,000 random points with Landsat 7 and Landsat 8 NDVI values 

obtained at the same locations. Next, linear regression was used to predict Landsat 7 NDVI with 

Landsat 8 NDVI (fig. B1; Area 2). There were significant differences between Landsat 7 and 8 

NDVI values before conversion, but the difference became very small after converting Landsat 8 

to Landsat 7 NDVI values with the following regression equation: 

NE_pre = -0.239 + (1.041 × NC) (B1) 

where 

 NE_pre is predicted Landsat 7 NDVI and 

NC is Landsat 8 NDVI.

This simple equation explained 94 percent of the variability between the two NDVI 

datasets and successfully converted Landsat 8 NDVI to Landsat 7 NDVI. The means and 

medians of the NDVI locations (n=5,000) were quite different before conversion, but almost 

identical following conversion (fig. B3A). The mean Landsat 8 NDVI was 0.2269 and the mean 

Landsat 7 NDVI was -0.0029 prior to conversion, but following conversion with equation B1, 

the Landsat 8 mean NDVI was -0.0028. 

Converting Landsat 7 NDVI to Landsat 5 NDVI

An almost identical procedure was used to convert between Landsat 7 and 5 NDVI values 

as those used to convert between Landsat 8 to 7 NDVI values except that Landsat 5 was not 

converted to TOA. Through experimental conversions of NDVI between sensors, I learned that 

the flycatcher satellite model simulations only produced accurate predictions of habitat when I 

did not adjust Landsat 5 to TOA before applying a correction factor. The reason for this is that 

the original satellite model of flycatcher habitat did not apply TOA correction (Hatten and 

Paradzick, 2003); thus, converting it retrospectively and then applying an NDVI correction 

between sensors resulted in unexpected results, such as a halo effect surrounding patches. 

Therefore, it is not advised to apply TOA correction to Landsat 5 imagery when applying the 

NDVI correction. 
85 

An area along the upper Gila River was used when NDVI values were converted between 

Landsat 7 and 5 (fig. B1; Area 1). After generating 5,000 random points in Area 1 and attributing 

them with Landsat 7 NDVI (corrected to TOA) and Landsat 5 (not adjusted to TOA) NDVI

values, a linear regression explained 0.93 percent of the variability in NDVI (fig. B3B) with the 

following equation: 

NT_pre = 0.102 + (0.897 × NE_pre) (B2) 

where 

 NT_pre is predicted Landsat 5 NDVI, and

NE_pre is predicted Landsat 7 NDVI. 

The means and medians of the two NDVI layers were quite different before conversion, 

but were almost identical after conversion (fig. B3B). Specifically, the mean NDVI value for 

Landsat 7 was -0.0858 compared to the Landsat 5 mean of 0.0247, but the Landsat 7 mean NDVI 

was 0.025 after conversion. 
86 

Figure B1. Map showing locations of areas 1 and 2 used in normalized difference vegetation index (NDVI) 

adjustments between sensors. Areas 1 and 2 were locations where 5,000 random points along the river 

channels were used to adjust Landsat 8 NDVI values to match Landsat 7 and 5 NDVI values. Project Area 

in this figure refers to a previous modeling effort where some of the methodologies presented in this report 

were originally developed (see Hatten, 2014).
87 

Figure B2. (A) Relation between Landsat 8 and 7 normalized difference vegetation index (NDVI), after 

adjusting to top-of-atmosphere [TOA] reflectance; and (B) relation between Landsat 7 and 5 NDVI, after 

adjusting to TOA reflectance. Landsat 8 and 7 images (A) were acquired in 2013, and Landsat 7 and 5 

scenes (B) were acquired in 2011. nc, Landsat 8 NDVI; ne, Landsat 7 NDVI; nt, Landsat 5 NDVI. NDVI 

values in (A) were acquired in Area 1 (fig. B1) while NDVI values in (B) were acquired in Area 2 (fig. B1). 

The numbers 3537 or 3637 refer to the different path/row numbers of each scene.

A

B
88 

Figure B3. (A) Distribution of Landsat 8 normalized difference vegetation index (NDVI) before 

(nc3637_131) and after (nc3637_131b) conversion to Landsat 7 (ne3637_131) NDVI, using equation B1; 

and (B) distribution of Landsat 7 NDVI before (ne3537_11174r) and after (nt3537_nt) conversion to 

Landsat 5 (nt3537_11166) NDVI, using equation B2. Landsat image ne3537_11174r was converted to topof-atmosphere (TOA) reflectance before NDVI calculations were made, whereas nt3537_11166 was not 

converted to TOA.

A

B
Publishing support provided by the U.S. Geological Survey

Science Publishing Network, Tacoma Publishing Service Center

For more information concerning the research in this report, contact the

 Director, Western Fisheries Research Center

U.S. Geological Survey

6505 NE 65th Street

Seattle, Washington 98115 

http://wfrc.usgs.gov/
ISSN 2331-1258 (online)

http://dx.doi.org/10.3133/ofr20161120

Hatten—A Model of Southwestern Willow Flycatcher Habitat and Effects of Tamarisk Leaf Beetles, Southwestern United States—Open-File Report 2016–1120