United States

Department of

Agriculture

Forest Service

Pacific Northwest

Research Station

General Technical 

Report

PNW-GTR-841 

August 2011

A Landscape Model for 

Predicting Potential Natural 

Vegetation of the Olympic 

Peninsula USA Using 

Boundary Equations 

and Newly Developed 

Environmental Variables

Jan A. Henderson, Robin D. Lesher, David H. Peter, 

and Chris D. Ringo
The Forest Service of the U.S. Department of Agriculture is dedicated to the principle of 

multiple use management of the Nation’s forest resources for sustained yields of wood, 

water, forage, wildlife, and recreation. Through forestry research, cooperation with the 

States and private forest owners, and management of the national forests and national 

grasslands, it strives—as directed by Congress—to provide increasingly greater service 

to a growing Nation.

The U.S. Department of Agriculture (USDA) prohibits discrimination in all its programs and 

activities on the basis of race, color, national origin, age, disability, and where applicable, sex, 

marital status, familial status, parental status, religion, sexual orientation, genetic information, 

political beliefs, reprisal, or because all or part of an individual’s income is derived from any 

public assistance program. (Not all prohibited bases apply to all programs.) Persons with 

disabilities who require alternative means for communication of program information (Braille, 

large print, audiotape, etc.) should contact USDA’s TARGET Center at (202) 720-2600 (voice 

and TDD). To file a complaint of discrimination, write USDA, Director, Office of Civil Rights, 

1400 Independence Avenue, SW, Washington, DC 20250-9410 or call (800) 795-3272 (voice) 

or (202) 720-6382 (TDD). USDA is an equal opportunity provider and employer.

Authors

Jan A. Henderson was area ecologist, Robin D. Lesher is ecologist, and Chris D. Ringo

was geographic information system analyst, Pacific Northwest Region, Mount BakerSnoqualmie National Forest, 2930 Wetmore Avenue, Suite 3A, Everett, WA 98201; David 

H. Peter is an ecologist, Pacific Northwest Research Station, Forestry Sciences Laboratory, 

3625 93rd Avenue SW, Olympia, WA 98512-9193. Henderson is currently located in 

Edmonds, WA. Ringo is currently located in Mount Vernon, WA.
Abstract

Henderson, Jan A.; Lesher, Robin D.; Peter, David H.; Ringo, Chris D. 2011.

A landscape model for predicting potential natural vegetation of the Olympic 

Peninsula USA using boundary equations and newly developed environmental 

variables. Gen. Tech. Rep. PNW-GTR-841. Portland, OR: U.S. Department of 

Agriculture, Forest Service, Pacific Northwest Research Station. 35 p.

A gradient-analysis-based model and grid-based map are presented that use the potential 

vegetation zone as the object of the model. Several new variables are presented that 

describe the environmental gradients of the landscape at different scales. Boundary 

algorithms are conceptualized, and then defined, that describe the environmental 

boundaries between vegetation zones on the Olympic Peninsula, Washington, USA. The 

model accurately predicted the vegetation zone for 76.4 percent of the 1,497 ecoplots used 

to build the model. Independent plot sets used to validate the model had an accuracy of 

82.1 percent on national forest land, and 71.8 percent on non-national-forest land. This 

study demonstrated that a model based on boundary algorithms can be an alternative 

to regression-based models for predicting landscape vegetation patterns. Potential 

applications of the model and use of the environmental gradients to address ecological 

questions and resource management issues are presented.

Keywords: Environmental gradient model, potential natural vegetation, vegetation zone 

map, Olympic Peninsula USA, boundary equations, environmental gradients. 
Summary

The concept of potential natural vegetation (PNV), the plant community reflecting the environmental capability of a land area, is a core concept in plant ecology and natural resource 

management. Rigorous, consistent, and validated potential vegetation mapping has been 

technically elusive, while at the same time it has remained a persistent need for land management agencies and a potentially valuable tool for stratifying and informing ecological 

studies. Gradients of predictive variables are sought to facilitate modeling and mapping of 

PNV. This paper presents a new way to model and map PNV to help meet these needs. We 

also present a new set of variables that describe the environmental gradients of a landscape 

at different scales. 

A gradient-analysis-based model and grid-based map are presented that use the potential vegetation zone as the object of the model and the underlying environmental variables 

as the primary subject of research. Boundary algorithms were conceptualized, and then 

defined, that described the environmental boundaries between vegetation zones on the 

Olympic Peninsula, Washington, USA. These boundary equations took the form of elevation of the boundary as the dependent variable, with components of the environment representing the independent variables. This was contrasted with the traditional regression-based 

approaches, which focus on central tendencies of the target ecological entity. 

New environmental variables developed for use in this model included “total annual 

precipitation at sea level,” “mean annual temperature at sea level,” “fog effect,” and “cold 

air drainage effect.” “Topographic moisture” and “temperature lapse rate” were newly 

calculated or represent new scaling for these concepts. Other more classic variables used 

as independent variables were “aspect” and “potential shortwave radiation.” 

The map produced by the study displayed the predicted pattern of eight potential 

vegetation zones. The model accurately predicted the vegetation zone for 76.4 percent of the 

1,497 ecoplots used to build the model. An independent validation set of 155 plots showed 

an accuracy of 77.4 percent.

This study demonstrated that a model based on boundary algorithms can be an alternative to regression-based models. In addition, a new set of environmental variables based 

on first- or second-order derivatives of weather station data or new hypotheses of landscape 

processes were used to successfully predict landscape patterns of potential vegetation. 

These maps of potential vegetation and the environmental variables can provide the 

foundation for subsequent scientific and ecological studies and for addressing resource 

management questions.
Contents

1 Introduction

2 Background, Concepts, and Definitions

2 Potential Natural Vegetation

3 The Plant Association and Vegetation Zone

5 Gradient Analysis

6 Mapping Considerations

7 Objectives

7 Study Area

8 Methods and Analysis

8 Field Methods

9 Analysis and Model Building

21 Model Validation

21 Results

21 New Environmental Variables

22 Model Building

23 Boundary Equations

24 Map of Potential Vegetation Zones for the Olympic Peninsula

27 Discussion

31 Acknowledgments

31 Acronyms

31 English Equivalents

31 References
GENERAL TECHNICAL REPORT PNW-GTR-841

ii

This page was intentionally left blank.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

1

Introduction

The relationships between vegetation and environment are 

the core of the field of plant ecology. Understanding these 

relationships is the subject of most papers in the field today. 

Some studies use these relationships to map the distribution 

patterns of target organisms or communities (Ohmann and 

Gregory 2002, Rehfeldt et al. 2006), and others display the 

relationships in statistical tables or abstract space (McKenzie et al. 2003, Rehfeldt et al. 2008). Most studies are limited 

to using existing and easily measured or derived environmental variables. 

This study used an innovative approach to predict vegetation zones (VZ) across a complex landscape, and several 

newly defined environmental variables are introduced. A 

hypothesis is presented that boundary equations can be 

used in lieu of more traditional regression-based algorithms 

(e.g., Ohmann and Gregory 2002, Rehfeldt et al. 2006). This 

hypothesis was tested by producing a map of the potential 

vegetation zones of the study area and by comparing it to an 

independent set of plots for validation. 

Potential vegetation is the object of this study because 

it is more closely related to the environment, as opposed to 

existing vegetation, which is strongly affected by seral stage 

and disturbance. Potential vegetation is the reference point 

or context for describing successional relationships and correlations between the vegetation and environment (Henderson et al. 1992). Potential vegetation is used in science and 

natural resource management for stratifying land relative 

to the environment and by informing questions regarding 

succession and growth potential.

The part of the environment that directly affects individual organisms has been called the “operational environment” (Mason and Langenheim 1957) as opposed to the full 

range of possible or measurable environmental variables. 

The operational environment is currently known only in a 

general sense; much of its apparent or correlative effects are 

thought to be related to some aspects of moisture, temperature, and light. 

Relative to terrestrial vegetation, moisture and temperature are highly complex forcing factors that vary in amount, 

phase, timing and quality in a wide array of possible 

combinations. These combinations of factors result in a 

unique environment and vegetation for nearly every place 

in the world. Part of the aim of plant ecology is to detect the 

effects of particular aspects of the moisture and temperature 

regimes, plus other factors including wind, light, and nutrients, while at the same time trying to uncover the complex 

ecological interactions, compensations, and variations over 

time and space.

The “operational” effects of this complex environment 

affect each individual organism independently (Mason and 

Langenheim 1957), in a corollary to the “individualistic” 

theory of vegetation organization (Gleason 1926). The sum 

of all such individual effects over time plus the effects of 

intra- and inter-specific interactions controls the composition of the terrestrial community. The integration of all 

such effects drives the structure and function of the ecosystem. Whereas the operational environment as a concept is 

relatively easy to identify in general terms, its parameters 

are poorly quantified for most species. However, the concept 

allows us to target for study those elements of the environment that appear to control or affect individual organisms, 

and by extension, communities and ecosystems. 

Understanding the operational environment and resulting distributions of species and communities has many 

practical applications. It can contribute to the mapping of 

species and communities, a necessary tool in the management of natural resources, biodiversity, and the conservation 

of biotic communities. 

Maps of the spatial patterns of units of vegetation 

(either existing or potential), as well as underlying environmental variables (such as different components of the 

moisture and temperature regimes), have a variety of uses 

for land management. These include analysis of biodiversity 

and planning for species or community protection; planning for restoration of habitats for at-risk species (such as 

the northern spotted owl [Strix occidentalis caurina] or 

marbled murrelet [Brachyramphus marmoratus]); planning 

for the protection of individual species, particularly those at 

risk; managing for beneficial or economic species; planning 

for timber management objectives such as yield, rotation 

length, or intermediate harvests; use as a benchmark or 
GENERAL TECHNICAL REPORT PNW-GTR-841

2

stratification for disturbance, climate change, and carbon 

sequestration studies and other monitoring applications; 

predicting the potential occurrences of wildfire and the rate 

of ecosystem recovery following perturbation; or describing 

the ecosystem departure from a sustainable range within an 

area. Maps of potential natural vegetation (PNV) are key 

to understanding patterns of existing vegetation as well as 

potentials for growth, species distributions, fire occurrence, 

and stand response to treatments. This is because of the 

strong link between units of PNV and the environment.

Historically, patterns of potential or existing vegetation 

have been depicted using intensive field surveys and handdrawn polygon maps (Dodwell and Rixon 1902, Eyre 1980, 

Kuchler 1964, Pfister et al. 1977). More recent attempts use 

satellite images and geographic information system (GIS)-

based digital data combined with weather station records to 

portray vegetation patterns in a pixel-based format (Cibula 

and Nyquist 1987, Ohmann and Gregory 2002). These digital mapping efforts use an algorithm to link environmental 

data to reference data of species presence (or sometimes 

absence), known plant community patterns, or field samples 

of plant community classes. This general modeling approach 

has its roots in the concept of “gradient modeling” (Kessell 

1979) and is sometimes known as “predictive vegetation 

mapping” (Franklin 1995).

The purpose of this study is to build a predictive, 

environmentally based, spatially explicit, gradient model of 

PNV using boundary equations to distinguish model units. 

In addition, and as necessary, the purpose is to develop 

new predictive environmental variables that help define the 

operational environment. 

Background, Concepts, and Definitions 

Potential Natural Vegetation

Tuxen is credited with the original definition of PNV, as 

“...the vegetation that would become established if all successional sequences were completed without major natural or direct human disturbances under present climatic, 

edaphic, and topographic condition” (adapted from Tuxen 

[1956], as translated and cited by Mueller-Dombois and 

Ellenberg [1974] and Kuchler [1964]). The more modern 

term “climax plant community” (CPC) is synonymous with 

Tuxen’s potential natural vegetation. However, the concepts 

embodied in the terms climax plant community, potential 

natural vegetation, association, and succession have evolved 

since 1956, especially regarding the link between vegetation and the environment, the nature of the climax, and the 

nature of vegetation itself. 

The concept of the CPC (climax association sensu 

Daubenmire 1952, 1976; Daubenmire and Daubenmire 

1968) is a theoretical classification unit representing the 

hypothetical end point of succession. As most plant communities (at least in the study area) are in some developmental/

successional stage moving toward this hypothetical climax, 

such a classification of PNV (sensu Tuxen 1956) or CPC 

(association sensu Daubenmire 1952) is an interpretation of 

the potential for successional development of the site. The 

existing communities of any area differ inherently in terms 

of their successional stage and their environment. Projecting such existing stands forward in time to a theoretical end 

point implies either projecting each individual community 

by itself or projecting a set of communities with similar 

environments. In either case, the climax community type 

concept implies a set of similar communities with similar 

environments that follows a repeatable or predictable 

successional pathway. 

A classification of the climax communities of an area 

(e.g., Henderson et al. 1989), without its environmental

context, means relatively little in terms of describing the 

community, its successional stages, and developmental 

potential or potential composition. We believe the link 

between the environment and the CPC has been implicitly 

acknowledged since the time of Daubenmire and Tuxen, 

and that it is appropriate to formally acknowledge the link 

between the environment and the nature of the CPC. The 

nature of the climax community is a function of the genetic 

diversity of the flora and its various and complicated interactions with the environment. It is continually changing over 

time as Earth’s climate and air and water circulation patterns change. Therefore, the nature of the CPC is dynamic 

through time. A classification or model at one point in time 

can only describe the potential vegetation for the climate of 

the time.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

3

The concept of potential natural vegetation used here 

refers to the set of communities, including all successional 

stages that can exist within the spatial and environmental 

bounds set by the extent of the CPC. The bounds of the CPC 

are controlled by the operational environment acting on 

the flora. The environmental bounds of a CPC type are the 

same as for all successional communities that may precede 

it, assuming the regional climate is stable. Therefore, PNV 

is identified or classified by the nature of the climax or 

projected climax community and is determined by the environment. This is not unlike some of the concepts embodied 

in the terms “habitat type” and “zone” (Daubenmire 1952, 

Daubenmire and Daubenmire 1968). 

The Plant Association and Vegetation Zone

The plant association (i.e., the CPC type as used here) has 

been used as the basic unit of the PNV classifications for 

the Western United States (e.g., Brockway et al. 1983; 

Daubenmire 1952; Daubenmire and Daubenmire 1968; 

Diaz et al. 1997; Franklin et al. 1988; Henderson et al. 1989, 

1992; Lillybridge et al. 1995; Topik et al. 1986; Williams 

and Lillybridge 1983). However, some authors (Daubenmire 

and Daubenmire 1968, Pfister and Arno 1980, Pfister et 

al. 1977) also used the concept “habitat-type” for the land 

area (mapping unit) where the plant association can develop 

(Daubenmire 1952, Daubenmire and Daubenmire 1968). 

Their emphasis on the land unit where a particular plant 

association occurs was used to emphasize the environment 

of the site. Plant associations can be aggregated into broader 

classification units of PNV called the vegetation series. 

Daubenmire’s [plant] association “is a concept embodying those characters of all actual stands among which differences in species composition are attributable to historic 

events or chance dissemination rather than to the inherent 

differences in environments” (Daubenmire 1952). In this 

somewhat convoluted definition, Daubenmire seems to try 

to establish the link not only between the climax community 

and the environment but to establish the current stand as 

part of the environmental potential of the site. It is not clear, 

however, if considerations of the breadth of environmental 

conditions within the range of an association, or the effects 

of factor compensation were included in his concept of 

(climax plant) association. 

Daubenmire’s “habitat type” linked the climax plant 

community (association) to the area of land where it can 

occur. By identifying the area of land where a single unit of 

the CPC type can occur, he explicitly linked the environment of that site to the developmental potential of the vegetation. Thus by extension, he appears to acknowledge the 

natural range in environments from one part of the habitat 

type to another. Other terms have been used to denote similar relationships such as the “range site” (Hironaka 1985), 

“site type” (Cajander 1926), or “ecological site” (Pellant et 

al. 2005, Smith et al. 1995). More recently, vegetation, soil, 

climate, and geology have been linked in a single ecological unit called the terrestrial ecological unit inventory 

(Winthers et al. 2005). Daubenmire (1952) also identified a 

broader unit of land he called the “zone,” which he defined 

as the area occupied or potentially occupied by a closely 

related group of (climax plant) associations. The National 

Vegetation Classification system also used the term plant 

association, but there it referred to existing vegetation rather 

than potential vegetation (Jennings et al. 2003). 

The VZ is the mapping unit of potential vegetation 

used in this paper. It approximates the area occupied by a 

vegetation series in scope and scale (as the terms are used 

in Daubenmire and Daubenmire 1968; Franklin and Dyrness 1973; Franklin et al. 1988; Henderson et al. 1989, 1992; 

Pfister et al. 1977) and is similar to the zone of Daubenmire 

(1952) and VZ of Franklin and Dyrness (1973). It is also 

similar in scope and ecological range to the “climax vegetation types” of Whittaker (1956) for the Great Smoky Mountains, and the 10 units of vegetation described and named by 

Whittaker and Niering (1965) for the Santa Catalina Mountains of Arizona. They variously described their 10 units of 

vegetation as vegetation types, zones, or coenoclines. The 

VZ is narrower in scope than the biotic communities of 

Brown et al. (1998) and Rehfeldt et al. (2006) or “potential 

vegetation type” of Ohmann et al. (2007).

The concept of the VZ used here is defined as “a unit 

of land where a single series or group of ecologically similar 

series predominates.” The characteristic or “zonal” series 
GENERAL TECHNICAL REPORT PNW-GTR-841

4

Figure 1—Location of the study area, the Olympic Peninsula, in western Washington state, USA. 

The area of the Olympic National Forest on the peninsula is shaded in darker gray.

(the taxonomic unit representing an aggregate of plant associations with the same climax indicator tree species) usually 

dominates within the broad geographic boundaries of a VZ, 

but azonal microsites may be occupied by vegetation of 

different physiognomic types and series. These sites usually 

represent the dry (e.g., a rock outcrop) or wet areas (e.g., a 

wetland) that could otherwise be seen as inclusions within 

the broad vegetation matrix characterized by a single series 

or small group of ecologically similar series. Subdivisions of 

the VZ are the plant association group (PAG) and the plant 

association (or habitat-type sensu Daubenmire).

An example of a VZ of ecologically similar but physiognomically different plant communities is the subalpine 

parkland zone. It is composed of tree-dominated series and 

several shrub- and herb-dominated series often in juxtaposition and forming a tightly integrated ecological system. 
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

5

Figure 2—Key to the eight potential vegetation zones on the Olympic Peninsula. This key works 

like a simple dichotomous key except that all the leads that would say “not as above” have been 

eliminated. Consider each pair of leads as in a normal dichotomous key, and read the lead as a question. If it is true, follow the lead to the right, if not, go to the next lead below it, and treat that line as 

a question, and so forth. 

 Land area above the upper limit of upright tree growth......................ALPINE ZONE

 Land area below alpine zone

Zonal vegetation nonforest and usually with scattered subalpine trees 

or tree islands...........................SUBALPINE PARKLAND ZONE

Zonal vegetation forest (with scattered nonforest communities in very 

wet or very dry habitats)

 Zonal potential vegetation characterized by ≥ 10 percent cover of Sitka 

spruce....................................................SITKA SPRUCE ZONE

 Zonal potential vegetation characterized by ≥ 10 percent cover 

mountain hemlock.......................MOUNTAIN HEMLOCK ZONE

 Zonal potential vegetation characterized by ≥ 10 percent cover 

Pacific silver fir..............................PACIFIC SILVER FIR ZONE

 Zonal potential vegetation characterized by ≥ 10 percent cover western 

hemlock and/or western redcedar..........WESTERN HEMLOCK ZONE

 Zonal potential vegetation characterized by ≥ 10 percent cover 

subalpine fir........................................SUBALPINE FIR ZONE

 Zonal potential vegetation characterized by ≥ 10 percent cover 

Douglas-fir............................................DOUGLAS-FIR ZONE

In a slightly different way, two or more very similar treedominated series can be aggregated into a single ecological 

VZ. An example of this is the western hemlock zone, which 

applies to areas where either the western hemlock series or 

the western redcedar series predominates.

There are eight VZs that occur on the Olympic Peninsula (fig. 1). These are identified in the key (fig. 2). These 

VZs are distinguished by the presence of key indicator tree 

species or certain community structural conditions related 

to dominant overstory species. This key is used to distinguish these eight VZs, and is based on a threshold of 10 

percent cover of the indicator tree species for the potential 

vegetation. The basis for the 10 percent threshold is given 

below. Plant names for the tree species follow the taxonomy 

of Hitchcock and Cronquist (1973). 

Gradient Analysis

The origin of gradient analysis of vegetation is from Whittaker (1956, 1967) and is based on the “individualistic” 

concept of Gleason (1926). Gradient analysis is the approach 

of relating gradients of plant populations and community 

characteristics to gradients of environment (Whittaker 1956: 

63). It seeks to understand the structure and variation of the 

vegetation of a landscape in terms of gradients of environmental factors, species populations, and characteristics of 

communities (Whittaker 1967: 207). Gradient modeling is 

the mathematical process of evaluating the links between 

environmental gradients on a landscape and their effects on 

communities and organisms (Kessell 1979). 

Various analysis techniques have been used in gradient analysis. Direct gradient analysis is usually applicable 
GENERAL TECHNICAL REPORT PNW-GTR-841

6

Figure 3—Distribution of key indicator tree species along the elevational gradient, with nodes 

separating different vegetation zones, from ecoplot data on the Olympic National Forest, 

Olympic Peninsula, Washington. Where “TSHE” is Tsuga heterophylla (western hemlock), 

“ABAM” is Abies amabilis (Pacific silver fir), “TSME” is Tsuga mertensiana (mountain hemlock), “WHZ” is western hemlock zone, “PSFZ” is Pacific silver fir zone, “MHZ” is mountain 

hemlock zone, and “PklZ” is subalpine parkland zone.

to patterns of individual species abundance along a simple 

or complex environmental gradient (ter Braak and Prentice 

1994, Whittaker 1967). Sometimes the direct gradients 

are complex or geographical gradients such as “elevation” 

(Dymond and Johnson 2002; Lookingbill and Urban 2005; 

Oksanen and Minchin 2002; Whittaker 1956, 1960).

In this study, direct gradient analysis is reserved for 

processes that use or generate ordinations of species abundance or frequency along known or quantified environmental gradients (fig. 3). Indirect gradient analysis is used to 

refer to any analytical method whose objective is to uncover 

or elucidate an unknown environmental gradient. Factor 

analysis refers to those analytical techniques used to determine which factors are best correlated between a known 

set of environmental factors and some quality of vegetation. 

Ordinations may be used in direct or indirect gradient 

analysis (ter Braak and Prentice 1994). All four of these 

approaches were used in this study. 

Mapping Considerations

One problem encountered while developing a pixel-based 

map of PNV is the resolution of scale. That is, transferring 

the scale of the field plots (fine) and classification (moderately fine) to a landscape (moderately coarse) via a model. 

Often the resolution of the types of vegetation on the ground 

is finer than the resolution of pixels (or polygons) used to 

portray them. Thus a single 90-m (0.81-ha) pixel can contain 

two or more fine-scale field plots of different community 

types or plant associations. 

Maps based on pixels are restricted to representing 

patterns of vegetation that are groups of pixels; therefore, 

the spatial pattern represented is necessarily coarser than 

the size of the individual pixel and usually much coarser 

than the scale of the field plots. The size of a pixel is often 

a function of the technology being used and the constraints 

of computer hardware and software used to represent them. 

Quantifying map accuracy between field plots of one scale 
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

7

and a map of pixels of another scale is difficult, both conceptually and practically. 

Another major challenge in developing a spatial model 

of PNV is having or developing appropriate mappable 

environmental variables that can be linked mathematically 

and conceptually to the dependent variable. Existing environmental variables have conventionally been those easily 

generated from weather stations or from digital elevation 

models (DEM) (Ohmann and Gregory 2002). Previous 

models have typically used aspect, slope steepness, elevation (as an independent variable), mean annual temperature 

(MAT), and total annual precipitation (in recent years) as 

quantitative variables, and slope position, soil, or geology 

have been used as qualitative variables. 

Any mapping exercise must necessarily generalize the 

detail of a landscape. Some practical or thought model is 

usually employed to classify the mappable elements of the 

landscape and to link the mapping algorithm back to the 

land. In practice, it is difficult to account for the finer scale 

of all inclusions and mosaics of vegetation in a complex and 

variable landscape, thus resulting in a model and a representative map that is more general than the actual landscape. 

Objectives

There were four objectives in this study: (1) create a map of 

VZs for the national forest land on the Olympic Peninsula, 

(2) create a map of VZs for the entire Olympic Peninsula, 

(3) find or develop a set of predictive and ecologically 

meaningful environmental variables, and (4) develop a 

predictive model of VZs using boundary equations and 

underlying environmental gradients. Development of a map 

of VZs of the Olympic National Forest required developing 

a spatially explicit model of VZs based on the available data 

while concurrently determining, finding, or developing the 

set of environmental variables to be used by the model. 

The model was then extrapolated to the entire Olympic 

Peninsula.

Study Area 

The area of study is the Olympic Peninsula in Washington 

state, USA (fig. 1). It is bounded on the west by the Pacific 

Ocean, on the north by the Strait of Juan de Fuca, on the east 

by Hood Canal, and on the south by an approximate line 

running southwest from the southern tip of Hood Canal, to 

Grays Harbor and the Pacific Ocean. The Olympic Mountains form much of the interior of the peninsula, centered 

about 123° 45 min. W. long. and 47° 45 min. N. lat. At the 

center of the Olympic Mountains lies the Olympic National 

Park (380 418 ha), which is surrounded on most sides by the 

Olympic National Forest (254 961 ha) (fig. 1). These two 

federal ownerships make up the bulk of the mountainous 

terrain of the peninsula. The lowlands are mostly held by 

state and private owners, plus six Indian reservations. Total 

land area of the study area is 1 490 773 ha. The elevation 

ranges from sea level to 2427 m at Mount Olympus. 

The Olympic Peninsula is known for its range of 

climates over a very short distance (Henderson et al. 1989, 

Peterson et al. 1997). Annual precipitation ranges from a 

total of 3353 mm at Quinault (western Olympics) to a low of 

432 mm at Sequim (northeastern peninsula), a distance of 

only 88 km. The western Olympic Range is characterized by 

the wettest and most maritime climate in the conterminous 

United States. The dry northeast corner is the driest place 

in western Washington and is an unusual mix of continental 

and dry-maritime climates.

The Olympic Range is geologically young and immature, originating near the end of the Pliocene Epoch. It consists of a core of sedimentary rocks, which are surrounded 

on three sides by an outer ring of mostly marine volcanic 

basalts (Tabor 1975). Mountain building is ongoing and is 

the result of the collision of the Juan De Fuca and North 

American plates and erosion, mostly owing to the numerous Pleistocene ice ages. Glacially carved U-shaped valleys 

characterize most of the Olympic Peninsula with some “V”-

shaped exceptions, especially in the south, indicating that 

recent glaciers were limited to the highest elevations and 

that water erosion has created much of the topography.

Vegetation is highly variable on the Olympic Peninsula 

owing to the strong maritime/rainshadow gradient from 

southwest to northeast and the rugged topography (Henderson et al. 1989). The Sitka spruce zone (SSZ) occurs at low 

elevations mostly along the Pacific coastline, where Picea

sitchensis (Bong.) Carr. (Sitka spruce) and Tsuga heterophylla (Raf.) Sarg. (western hemlock) dominate. Elsewhere, 
GENERAL TECHNICAL REPORT PNW-GTR-841

8

the western hemlock zone (WHZ) occurs at lower elevations 

where Pseudotsuga menziesii (Mirbel) Franco (Douglasfir) and western hemlock dominate. Thuja plicata Donn 

(western redcedar) may occur as an occasional species or 

codominant at low to mid elevations. The Pacific silver fir 

zone (PSFZ) occupies mid elevations, where western hemlock and Abies amabilis (Dougl.) Forbes (Pacific silver fir) 

usually codominate, along with small amounts of western 

redcedar and, occasionally, Douglas-fir as a long-persistent 

seral species. The mountain hemlock zone (MHZ) occurs 

at higher elevations, where Pacific silver fir and Tsuga 

mertensiana (Bong.) Carr. (mountain hemlock) codominate. Chamaecyparis nootkatensis (D. Don) Spach (Alaska 

yellow-cedar) occurs infrequently at mid to high elevations. 

The subalpine fir zone (SAFZ) occurs at higher elevations 

in the drier areas of the northeast Olympics where Abies 

lasiocarpa (Hook.) Nutt. (subalpine fir) dominates. The 

Douglas-fir zone (DFZ) also occurs to a limited extent in the 

northeast Olympics, at low to mid elevations and on southerly and drier topographic positions. Here Douglas-fir is the 

dominant climax species, and the environment is too dry for 

western hemlock or western redcedar. At upper elevations, 

the forest zones give way to the subalpine parkland zone 

(PKIZ) composed of tree islands and subalpine meadows 

at about 1067 m in the west and at elevations over 1524 m 

in the northeast part of the peninsula. Small areas of alpine 

zone occur at the highest elevations. 

Methods and Analysis

Field Methods 

Four types of field data were used to develop and validate 

the VZ model presented here. The first type is the “ecoplot” 

data collected as part of the U.S. Forest Service (USFS) 

Region 6 Ecology Program (Henderson et al. 1989, Henderson et al., unpublished data). Data collected on these fixedarea plots included, but were not limited to (1) location, (2) 

elevation, (3) direction of the sloping land (i.e., aspect), (4) 

steepness of the sloping land (slope), (5) topographic moisture (Henderson et al. 1989, 1992; Whittaker 1956, 1960; 

Whittaker and Niering 1965), (6) identity of each vascular 

plant species according to Hitchcock and Cronquist (1973) 

plus a measure of their individual cover, (7) plant association 

(Henderson et al. 1989), vegetation series, and vegetation 

zone. A total of 1,497 ecoplots were used in this study; these 

plots were well distributed in different watersheds and represent a wide range of ecological conditions on the Olympic 

National Forest. The locations of many of these plots were 

randomized near the center of accessible General Land Office sections (square miles) of land. This was the primary 

data set used to build the VZ model.

The second type of plot data used was a set of 447 check 

plots collected along road transects throughout the Olympic 

Peninsula. These plots occurred on all types of ownerships 

where access was available or where permission to enter was 

granted by the owner or administrator. Data on these plots 

included (1) location by Rockwell1 “PLGR” global positioning system (GPS), (2) elevation, and (3) VZ. Other data were 

extracted by intersecting the GPS location (if local error was 

less than 10 m) with GIS grids of topographic or environmental variables. These plots were used in model building 

to extend the ground-truthed part of the map to areas not 

previously sampled, including extensive areas of private or 

state-owned lands.

A third type of plot data used for model building was 

field plots (n = 1,051) contributed by the Olympic National 

Park (ONP). These data are represented by three sets called 

“ONP-Elwha” (n = 63) (A. Woodward, U.S. Geological 

Survey [USGS], unpublished data), “ONP Animals” (n 

= 431) (K. Jenkins, USGS, unpublished data), and “Hoh/

Dosewallips 1979” (n = 557) (J.A. Henderson and B.G. 

Smith, ONP, unpublished data). These field plots were taken 

using various field protocols from 1978 to 2002. Plot locations were extracted from field data, maps, or aerial photos. 

Vegetation zone identifications were done by National Park 

Service personnel or contractors and checked by the authors. 

The protocol used for the “Hoh/Dosewallips 1979” data set 

was similar to the field protocol used in the USFS ecoplots 

described above.

1

The use of trade or firm names is for reader information and does not imply endorsement by the U.S. Department of Agriculture of 

any product or service.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

9

A fourth set of independent plots was used to validate 

the final model. These plots were a subset of USFS Current 

Vegetation Survey (CVS) plots (n = 84) (USDA 1998) and 

Forest Inventory and Analysis (FIA) plots (n = 71) (Campbell et al. 2010). These systematically distributed permanent 

plots are maintained across federal, state, and private ownerships for the purposes of periodic national timber inventory 

mandated by Congress in the Forest and Rangeland Renewable Resources Planning Act (RPA) of 1974. These plots 

were first screened for accuracy of GPS location and then 

assigned to a vegetation zone. Those plots with field collected data sufficient for reliable placement into a VZ were then 

used as an independent validation sample to determine the 

accuracy of the final version of the VZ map of the Olympic 

Peninsula.

To summarize, the set of USFS ecoplots (n = 1,497) 

was used for model building and environmental variable 

identification. The check plots (n = 447) and ONP plots (n = 

1,051) were used to extrapolate the model to non-nationalforest lands. Two independent data sets were used for model 

validation, the CVS plots on national forest land (n = 84), 

and FIA plots on non-national-forest lands (n = 71). The total 

number of field plots used in this study was 3,150.

Analysis and Model Building

Overview of model building—

This study was done in three parts. Each part represented 

a different approach with different objectives. Part one set 

the broad environmental gradient context and identified the 

boundary nodes between the VZs (figs. 2 and 3). Part two 

was the development of the new environmental variables 

and adaptation of existing variables used in the model. 

Part three was the development of the boundary equations 

between adjacent pairs of VZs. 

Part one consisted of extensive field survey and classification of climax community types (associations sensu 

Daubenmire) (Henderson et al. 1989). It also used a direct 

gradient analysis of the tree component of the climax 

vegetation to help define boundaries between VZs of the 

Olympic National Forest (fig. 3), and was later expanded to 

the Olympic Peninsula (fig. 1). In these initial stages, tree 

cover measured on the 1,497 ecoplots was summarized by 

elevation classes from sea level to the highest peaks (2126 

m). This broad averaging of tree cover (fig. 3) over such a 

wide range of moisture and temperature gradients necessarily obscured some of the range in environmental conditions 

across the Olympic National Forest. Similar graphs were 

also done for different geographic areas, and for different 

precipitation and temperature regimes, here defined. Each 

revealed different patterns of tree cover distribution and 

elevation. However, for presentation, only the general form 

is given, as it represents the conceptual basis for determining the boundaries of the VZs.

The nodes identified along this complex gradient were 

used to identify the boundaries between different VZs in 

the study area (fig. 3). For most boundaries, potential presence of 10 percent cover of the indicator tree species at its 

warmer or drier boundary indicated the lower elevational 

boundary of the zone (Henderson et al. 1989). Exceptions to 

this are that vegetation physiognomy was used to distinguish 

between the upper forest boundary and the subalpine parkland zone, and tree form was used to distinguish between 

the subalpine parkland and alpine zones (fig. 2). Ten percent cover is understood to be an ecologically meaningful 

abundance threshold for zonal sites in the landscape. Below 

10 percent cover, the distribution of indicator tree species is 

sporadic and they tend to occur in azonal microsites. These 

nodes fixed the boundary between any two vegetation zones 

in real space as well as in environmental hyperspace, for any 

given combination of model variables. 

Once these nodes were established, they were used to 

assign VZ codes to each sample plot in the study. Because 

this required interpretation of the stand age and structure to determine whether at least 10 percent cover of the 

indicator tree species was present or potentially present in 

late-successional development, some knowledge of the tree 

species and the environments of the area was needed. For 

forest communities older than 300 years, the actual covers 

of late-successional tree species were used to represent the 

climax condition. Most plots were easily placed into one of 

the vegetation zones. A few plots were clearly transitional, 

and were assigned on best judgment, but some error was 

probably created at this time. A few plots, often because of 
GENERAL TECHNICAL REPORT PNW-GTR-841

10

their young age, could not be reliably assigned to a VZ, so 

they were omitted from further analysis. Plots with stand 

ages less than 10 years were not used.

Part two was the development of a set of environmental 

variables used to predict the VZ boundaries. Once the first 

variable was identified, the examination of residual errors 

led to the discovery of better or additional environmental 

variables. This was a complex process, but was initiated 

when the use of existing environmental variables and 

methods did not provide the means to adequately describe 

the pattern of VZs on the landscape or give the ecological 

understanding of the underlying causes of the distribution of 

vegetation. These variables included several environmental 

indices defined later in this paper including precipitation 

at sea level (PSL), mean annual temperature at sea level 

(MATSL), fog effect, topographic moisture (TM), and cold 

air drainage (CAD) effect. At that time, data sets to serve 

this purpose were not available at a resolution adequate for 

the study objectives (e.g., Daly et al. 2002, Hungerford 

et al. 1989).

Part three was the development of the boundary equations. We used graphical analysis of sets of plots from 

adjacent pairs of VZs to estimate the mathematical boundary between each such pair of VZs (fig. 4). The boundary 

line between the upper boundary of one VZ and the lower 

boundary of another was determined by the fit of an equation to achieve approximately equal numbers of errors above 

and below the line with a minimum of total errors. 

This iterative process of analysis was repeated at each 

step in the uncovering of each new environmental variable. 

Thus part two and part three were used iteratively to evaluate and develop each variable. 

These three parts culminated in a map of vegetation 

zones. Converting the mathematical equations to a spatial 

model began in the early 1980s by manually mapping the 

VZs in the field using topographic maps and field data. With 

the development of technology over subsequent years, this 

mapping process developed into a complex GIS computer 

application (the PNV Model) programmed in Arc macro 

language and processed in the Arc/Info Grid program (ESRI 

1982–2006). 

Definition and development of the model variables—

Precipitation at sea level—

Precipitation at sea level is a new variable presented here. 

It represents a precipitation regime for the landscape and is 

defined as the relative precipitation for an area with the 

direct effects of elevation removed, and is portrayed at the 

base elevation of “sea level.” Precipitation varies greatly 

over the Olympic Peninsula and is not the same at a given 

elevation, as was observed when contrasting the southwest 

and the northeast parts of the peninsula. However, where 

PSL is the same, precipitation will be the same at a particular elevation throughout that area. Precipitation at sea level 

values describe proportionately how much more or less 

precipitation can be expected across a landscape or geographical area. Thus the pattern of total precipitation can be 

described as two vectors: (1) the vertical, related to elevation 

and (2) the horizontal described by the pattern of PSL. It 

is, therefore possible to calculate actual precipitation at any 

point on the landscape by knowing its PSL and elevation as 

will be shown later. 

Precipitation at sea level can be regarded as an index of 

relative precipitation potential, similar to how forest stand 

growth is often indexed by the projected height of dominant 

and codominant trees at a base age, i.e., the “site index.” 

Developing an index that allowed the spatial representation of the relationship of precipitation to a base elevation 

(i.e., sea level) was a breakthrough that greatly facilitated 

the modeling of potential vegetation over the landscape. 

Because the model identified elevations of the boundary 

between VZs, it was critical to identify areas where the 

relationship of precipitation to elevation remained constant. 

The development of PSL involved a series of successive 

iterations, and provided insight into both the necessity for 

and the utility of this variable. 

The history of the development of the PSL variable 

began in 1982 with the initial attempt to field map the VZs. 

This process began with analysis to map the boundary 

between the western hemlock zone and the Pacific silver fir 

zone. These two zones were represented by the greatest 

amount of data and covered the greatest area across the 

Olympic National Forest. 
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

11

Charts were first developed to depict the aspect/elevation patterns of the lower boundary Pacific silver fir zone 

(LBPSFZ) from one geographical area to another across the 

Olympic National Forest. For this effort, soil temperature 

measurements recorded on USFS ecoplots were used to 

estimate the magnitude of elevation effect from northerly 

to southerly aspects. This analysis showed a change from 

north to south aspect at the lower boundary of the PSFZ 

equivalent to a change of about 200 m elevation (Henderson 

et al. 1989). Also, it indicated that the greatest effect was at 

aspects of about 30 degrees (where the LBPSFZ was about 

100 m lower) and at 210 degrees (where the LBPSFZ was 

about 100 m higher elevation). Initially, sample plots from 

a local area in the South Fork Skokomish River watershed 

(southeast Olympics) were plotted on a simple graph with 

aspect on the x axis and elevation as the y axis (fig. 5A). A 

sine function was fit to represent the lower boundary of the 

PSFZ and the upper boundary of the WHZ. A similar curve 

was fit to the PSFZ/MHZ boundary.

The relationship represented by the curves in this graph 

(fig. 5A) was applied to the South Fork Skokomish watershed 

by hand-drawing the lower boundary of the PSFZ on USGS 

topographic maps. It correlated well with the vegetation pattern observed in the field. However, when these curves were 

used to make a similar map for the Wynoochee watershed 

to the west (toward the Pacific Ocean), the relationship was 

weak. In this case, the predicted LBPSFZ (curves from fig. 

5A) was too high in elevation to fit the data (fig. 5B). The 

residual errors between the predicted function and the actual 

elevation of the LBPSFZ were computed and the boundary 

curves were recalculated. Similarly, when nearby watersheds (North Fork Skokomish and Hamma Hamma) were 

examined to the northeast, the original function (fig. 5A) 

predicted a boundary that was too low compared to the data 

in fig. 5C. Figures 5B and 5C show the new lower boundary lines for PSFZ and MHZ that fit the data for these two 

geographic areas. Note the difference between the boundary 

curves from fig. 5A compared to the boundary curves in 

fig. 5B (wetter area) and fig. 5C (drier area).

This trend was noted across the entire Olympic Peninsula. It was not known at the time which environmental 

factor(s) was causing the aspect/elevation shift. This did not 

Figure 4—Vegetation zone plots displayed by elevation and precipitation at sea level, with 

boundary equation lines for lower boundary Pacific silver fir zone and lower boundary 

mountain hemlock zone. Where “MHZ” is mountain hemlock zone, “PSFZ” is Pacific 

silver fir zone, “WHZ” is western hemlock zone, “LBMHZ” is lower boundary mountain 

hemlock zone, and “LBPSFZ” is lower boundary Pacific silver fir zone.
GENERAL TECHNICAL REPORT PNW-GTR-841

12

Figure 5A

Figure 5B

Figure 5—Vegetation zone plots displayed for three small 

geographic areas in the southern Olympic Mountains, 

with boundary curves between vegetation zones. 

“MHZ”is mountain hemlock zone, “PSFZ” is Pacific 

silver fir zone, “WHZ” is western hemlock zone, 

“LBMHZ” is lower boundary mountain hemlock zone, and 

“LBPSFZ” is lower boundary Pacific silver fir zone. Area of 

graph 5A (South Fork Skokomish watershed) has average 

precipitation at sea level (PSL) of 2413 mm and mean 

annual temperature at sea level (MATSL) of 11.0 °C. Graph 

5B is for an area west (Wynoochee watershed) of graph 5A, 

with average PSL 2896 mm and MATSL 10.8 °C. Graph 

5C is for an area northeast (North Fork Skokomish and 

Hamma Hamma watersheds) of graph 5A, with average 

PSL 1753 mm, and MATSL 11.2 °C.

Figure 5C
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

13

prevent the empirical and systematic field mapping of the 

LBPSFZ across the Olympic National Forest. This serpentine boundary line traversed the area, up and down in 

elevation in response to many, as yet unknown environmental factors. It was lower in areas with wetter soils, higher 

precipitation, north aspects, concave slopes, colder areas, 

or in valleys with cold air accumulations. 

Areas where the LBPSFZ was located at the same 

arbitrary elevation were defined as “environmental zones” 

(Henderson et al. 1989: 40, fig. 24). These zones represented 

geographical areas where the environment was judged to be 

similar and were originally identified by codes from zero to 

12 where each zone is represented by the increase in elevation of the LBPSFZ of about 80 m. These geographical areas 

were used to analyze the upper boundary of the PSFZ and 

to define curves for it. Similarly upper and lower boundary 

“aspect/elevation” functions were defined for each boundary 

in the area. This tool was used to make the first preliminary 

map of vegetation zones for the area (Henderson et al. 1989: 

45, fig. 25).

The correlation between the “environmental zones”

and the existing precipitation map for the peninsula 

(Henderson et al. 1989: 35, fig. 18) was obvious. However, 

evaluation of apparent differences along the coast and other 

discrepancies inland suggested that other variables were 

also driving these “environmental zones” besides total annual precipitation. Some evidence pointed to an interaction 

between precipitation and elevation, other evidence pointed 

to an interaction with fog or some other maritime factor.

A number of transformations of the precipitation 

data were attempted to fit a better function between VZ 

boundaries and environmental zones. The relationship was 

somewhat obscured, partly because of the small number of 

weather stations on the peninsula and their aggregation at 

low elevation near saltwater. It was not until a similar study 

was done in the North Cascade mountains of Washington 

that the relationship between precipitation and elevation (to 

be called “precipitation at sea level”) was formalized.

The new variable “precipitation at sea level” (PSL) 

was generated by graphing the measured annual precipitation at weather stations against the elevation. At first, a 

single linear regression was fit to each of a number of sets 

of nearby weather stations across western Washington. 

Linear functions were used initially until it became apparent that a quadratic function would better fit the weather 

station data of precipitation relative to elevation. As there 

were few stations at higher elevations, this trend could not 

be established with any certainty. However, it was assumed 

that the rate of precipitation increase gradually decreased 

with elevation approaching an asymptote at the elevation of 

the highest ridges, and then declining slightly as air thinned 

with increasing altitude, as represented in equation 1. From 

this analysis, a single equation was generated to project the 

measured annual precipitation at each weather station back 

to an elevation equal to zero (eqn. 1): 

TAP = PSL + (A + B × PSL) × ELEV + 

(C + D × PSL) × ELEV 2

 (1)

Where

TAP = total annual precipitation (mm)

PSL = precipitation at sea level (mm) 

A = 3.07916667E-02

B = 5.18339895E-04

C = - 5.17648321E-06

D = -7.19779665E-08

ELEV = elevation in meters

Although PSL was originally solved graphically with 

a set of curves using equation 1 with different intercept 

values, it can be solved directly using equation 2. Equation 

1 was rewritten as follows (eqn. 2) to predict PSL for each 

weather station: 

 PSL = 1 / {M+N / TAP + [1 / (P+Q × TAP)] × ELEV + 

 [1 / (R+S × TAP)] × ELEV 2} (2)

Where

PSL = precipitation at sea level (mm)

M = 3.14066E-07

N = 0.99908387

P = -139015.01

Q = 1922.4891

R = 561677190

S = -13825146

TAP = total annual precipitation (mm)

ELEV = elevation in meters
GENERAL TECHNICAL REPORT PNW-GTR-841

14

These PSL values for each weather station on the 

Olympic Peninsula were calculated and plotted on a map 

and interpolated to generate an isoline map of PSL (fig. 6). 

Precipitation at sea level values for each pixel in the study 

area were derived by triangulation on these isolines. This 

map of PSL was later substituted for “environmental zones” 

and became the basis for the VZ model. However, the 

original pattern of environmental zones was retained as 

the starting point for drawing the PSL lines. Later as this 

analysis was extended throughout the states of Washington 

and Oregon, the PSL isolines were refined. 

The PSL isolines were then calibrated so that the model 

correctly predicted the average total annual precipitation 

for each station to less than 1.0 percent error and the sum of 

deviations was zero. Later a topographic/geographic analysis 

was done to determine on-shore flow and rain-shadow patterns. This was used to refine the original coarse pattern 

of PSL lines to better reflect the effects of topography and 

Figure 6—Precipitation at sea level (PSL) for the Olympic Peninsula. Contour lines represent 

PSL in millimeters. PSL values for each pixel in the study area were derived by triangulation 

on these curves.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

15

geography on the flow of clouds, storms, fog, and rain across 

and around the Olympic Peninsula (fig. 6). 

Once the environmental variable of PSL was developed, 

and the base model with elevation of the VZ boundary as the 

dependent variable and PSL as the independent variable was 

derived, the next step was to evaluate the residual errors of 

VZ prediction both mathematically and spatially to generate 

patterns of deviation. This was followed by hypothesizing 

and exploring possible environmental/ecological relationships in the landscape that could account for these errors. 

The successful development and application of the PSL 

concept led to the exploration of the idea that other new 

environmental variables related to the patterns of distribution of PNV could be developed and incorporated into the 

PNV Model application.

The next step in the process started by evaluating a 

relationship borrowed from Whittaker (1956) called “topographic moisture,” and its effect on the distribution of the 

PSFZ and later all VZs. Eventually four new variables in 

addition to PSL (topographic moisture, fog effect, MATSL 

cold air drainage effect), plus an adaptation of potential 

shortwave solar radiation (Kumar et al. 1997) were evaluated and added to the model.

Topographic moisture—

This is a quantification of the “topographic moisture gradient” concept of Whittaker (1956, 1967, Whittaker and 

Neiring 1965). Whittaker (1967: 218) presented species data 

along a scale from 1 to 13 representing “mesic to xeric” 

conditions, but without otherwise quantifying the gradient. Whittaker noted that the units of the scale were groups 

determined by ordination.

Conceptually, the topographic moisture (TM) variable 

as used here is this: as soil water moves downslope under 

the influence of gravity, it tends to move faster on steeper 

and more convex slopes and slow down or accumulate on 

gentle and concave slopes. This TM index is intended to 

scale the relative hillslope drainage (e.g., Beven and Kirkby 

1979). Therefore in areas with higher TM values, there is 

potentially more soil water availability (soil texture and 

other impediments being otherwise equal), and in areas with 

lower TM values, there is the topographic potential for less 

available water. In contrast to Whittaker’s (1956, 1967) use, 

TM as used here does not include effects of soil texture or 

aspect. Also it does not take into account potential subsurface waterflow at bedrock surface or at textural discontinuities in the soil or seeps or springs, as these features are 

mostly unmapped in the study area. 

For the present study, a scale of values from 1 to 9 was 

developed to represent the range of topographic moisture 

potentials for the landscape. The values ranged from 1 

(absolute driest topographic condition) to 9 (absolute wettest 

topographic position, i.e., open body of water), where the 

value of 5.0 represents the modal topographic moisture 

value for the landscape. These values were assigned to 

sample plots in the field (ecoplots) (Henderson et al. 1989) 

to use in calibrating the effect of this variable on the lower 

boundary of the PSFZ in western Washington. 

A landscape model was constructed to predict the 

topographic moisture value for pixels in the landscape. This 

algorithm used the “CURVATURE” function in the Grid 

Module of Arc Info Workstation (ESRI 1982–2006) plus 

adjustments for steepness of slope. These values were averaged to windows of 270 m and 900 m radius, to compute 

TM for the Olympic Peninsula. These two windows were set 

to reflect meso-topography at the local ridge scale and the 

broader scale of major ridges. 

Fog effect— 

The nature of the vegetation along the western peninsula 

indicated that the vegetation appeared to be responding to 

a wetter environment than was represented by precipitation alone (fig. 6) (Henderson et al. 1989: 35, fig. 18). Davis 

(1966) had suggested a relationship between spruce forests 

and coastal fog effects in Maine. A new effect, resembling 

the fog patterns presented by Stone (1936) or Vogelmann 

et al. (1968) was developed to reflect this “fog” effect on 

vegetation. Weather satellite images and field reconnaissance further indicated that some maritime effect appeared 

to contribute a surrogate of precipitation or to enhance the 

precipitation effect in this area. Isolines of fog effect along 

the coast of Washington and Oregon were drawn to reflect 

both cloudiness and the presumed effect of fog and dew on 

vegetation (fig. 7). Later, through direct observations and 
GENERAL TECHNICAL REPORT PNW-GTR-841

16

identification of topographic barriers, this map was refined 

further to reflect actual fog and low cloud penetration onto 

the land. This fog effect was weighted to reflect the average 

value of effective precipitation added to the site by fog drip; 

a number suggested by several studies in coastal California 

and Oregon. For this study we used a weighting of 508 mm 

of effective precipitation for one unit of fog effect. Azevedo 

and Morgan (1974) found up to 425 mm of summer fog drip 

in coastal redwood (Sequoia sempervirens (D. Don) Endl.) 

forests in an area with an undisclosed amount of summer 

precipitation. Oberlander (1956) measured 1524 mm of fog 

drip under an exposed tanoak (Lithocarpus densiflorus

(Hook. & Arn.) Rehder), and Dawson (1998) measured 447 

mm of annual fog drip in coastal redwood forests with annual precipitation of 1315 mm; fog drip amounted to 34 percent 

of the annual precipitation. Harr (1982) found 914 mm of 

fog drip in the Bull Run Watershed in the northern Oregon 

Cascades. Values of fog effect ranged up to 1524 mm in 

Figure 7—Fog effect for the Olympic Peninsula. Contour lines represent fog effect in millimeters equivalent of precipitation. Fog effect values for each pixel in the study area were derived 

by triangulation on these curves.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

17

some areas on the Pacific Ocean side of the peninsula to less 

than zero in the extreme rain-shadow areas of the northeast 

peninsula. Although this effect was correlated with fog, the 

ecological effects are perceived to be more complex, involving dew, condensation on leaves, stem flow, tree drip, etc.

Mean annual temperature at sea level— 

Mean annual temperature at sea level was constructed using 

the same approach and similar concept as for PSL. It is a 

new variable presented here that represents the temperature 

regime for a landscape, and is presented as an index of MAT 

at a base elevation of sea level (fig. 8). Linear regressions 

were done on weather station MAT and elevation using data 

from sets of proximal stations across northwestern Washington. The slope of these lines is the “temperature lapse 

rate” variable (fig. 9). The slopes of these lines were used 

to project the MAT of each weather station down to a base 

elevation of sea level. A digital data set of MATSL (fig. 8) 

Figure 8—Mean annual temperature at sea level (MATSL) for the Olympic Peninsula. Contour 

lines represent MATSL in degrees C. The MATSL values for each pixel in the study area were 

derived by triangulation on these curves.
GENERAL TECHNICAL REPORT PNW-GTR-841

18

was then generated using similar procedures as for PSL. 

Later the map of MATSL was calibrated using data from 

eight remote weather stations for the Olympic National 

Forest. 

Cold air drainage effect—

The cold air drainage (CAD) effect model used an algorithm 

to calculate a MAT value for each pixel in the landscape, 

using MATSL, elevation, lapse rate, and aspect. The aspect 

effect was eventually adjusted by adding the effects of shortwave (SW) solar radiation to the model. Conceptually, cold 

air at night, in winter, and in still air accumulates in concave 

areas represented in the digital terrain model. A proxy for 

air density was developed using modeled MAT in a cold-air 

flow model. This proxy variable was called “temperature 

units” and was calculated as the maximum MAT for the 

study area plus 1.1 °C less the calculated MAT of each pixel. 

Figure 9—Temperature lapse rate for the Olympic Peninsula. Contour lines represent the rate at 

which temperature declines per 1000 m elevation in degrees C. Lapse rate values for each pixel 

in the study area were derived by triangulation on these curves.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

19

ArcInfo’s “FLOWACCUMULATION” function (ESRI 

1982–2006) was used to identify potential cold air drainage channels across the landscape. The flow accumulation 

module in ArcInfo Workstation was used to determine the 

flow direction and flow accumulation following the natural 

flow patterns across a landscape. This process is described 

by Jenson and Domingue (1988) and Chung et al. (2006). 

It results in the sum of the number of pixels that are calculated to be above the target pixel and to be within its direct 

flow pathway as determined by relative elevation and aspect. 

By repeating this for each pixel in the landscape, the model 

calculates the “flow-accumulated” sum of pixels for each 

pixel in the landscape. This resulting grid contained higher 

values where there was greater potential for cold air accumulation. Therefore the value being accumulated through 

the “FLOWACCUMULATION” function was related to 

the density of air and the accumulation of high CAD effect 

values in areas in cold air drainages and cold air lakes. Thus 

the CAD model reflected the downward convective movement and accumulation of denser and colder air.

The flow accumulation was then modified by a diffusion function that applied a series of “FOCALMEAN” 

(ESRI 1982–2006) functions to “spread” the accumulated 

“temperature units” across increasing circles of influence 

dependent on the topography. The more open and flat a terrain, the more the flow accumulation is spread out among 

adjacent pixels. This final grid of accumulated temperature 

units (ATU) was then compared to the measured change 

in elevation of selected VZ indicator species in a sample of 

locations across the state of Washington. A calibration was 

made for the elevational change of indicator species and 

areas of high cold air accumulation where elevation shift 

owing to CAD is a function of ATU. The downward shift 

in VZ indicator species was calibrated to the ATU and a set 

of curves was developed where ATU and PSL were independent variables and elevation change expressed in meters 

for each VZ was the dependent variable. After repeated 

experiments in calibrating this function to the landscape, a 

final model was achieved. The output of the cold air drainage effect model was the estimated lowering of temperature 

in CAD accumulation areas expressed in terms of units of 

elevation effect on vegetation (fig. 10). A further calibration was made to calculate absolute temperature change (in 

degrees Celsius) relative to this predicted elevation change 

effect. 

Shortwave solar radiation—

The effects of potential radiation load for each pixel in a 

landscape owing to the angle of incidence and topographic 

blocking of solar radiation was used to supplement the role 

of aspect as a variable in the model. Although aspect as a 

predictor was adequate for mid-slope VZs and in mountainous terrain, for areas on ridgetops, in valley bottoms, and in 

broad gentle terrain it did not perform well. 

The solar radiation model of Kumar et al. (1997) was 

used to predict radiation for the area and was incorporated 

into the VZ model in combination with aspect. The final version of the VZ model includes a function for “aspect effect” 

that is: 

aspect effect = P × f (aspect) + (1 - P) × [(SW - 14930) / 

3380.4]

Where 

shortwave (SW) solar radiation in watts per square meter 

is transformed to units of aspect by the subfunction (SW - 

14930) / 3380.4 and where P and 1 - P are the proportions 

of the total aspect effect distributed between the variables 

aspect and SW. Thus a P of 0.5 would weight the direct 

effect of aspect equally with the transformed effect of solar 

radiation. This function is a complex variable that combines 

the effect of solar radiation and aspect by transforming solar 

radiation to units of aspect effect and weighting the two 

variables by P and 1 - P.

These environmental variables were originally added to 

the VZ model in the order they were developed or discovered, but the final step of model development evaluated their 

relative contribution. Thus variables were added to the final 

model in order of their effect on accuracy of VZ prediction 

of the field plots used to build the model.

Vegetation zone model building—

The model presented here used a different approach than 

other attempts to develop a landscape model of PNV. This 

model was based on boundary algorithms that predict the 
GENERAL TECHNICAL REPORT PNW-GTR-841

20

lower and upper elevational limits of each VZ, approximating an ecological envelope approach, rather than estimating 

the likelihood (e.g., Ohmann et al. 2007, Robertson et al. 

2001) of each pixel belonging to each potential class of PNV. 

The VZ model began with the lower boundary of the PSFZ 

equation suggested in figures 4 and 5 and is given as equation 3: 

 VZelev = A + Bx1 + Cx2 + Dx3 + ... Nxn (3)

Where 

VZelev = the elevation representing the node for the lower or 

upper elevational boundary of any VZ in the study area, A is 

the intercept term, B...N is the slope function coefficient or 

amplitude of the environmental variable x1...xn. In the case 

of aspect, the variable x is the transformed relationship of 

the sine of the aspect minus 120°. Furthermore, the dependent variable is the point in geographic space represented 

Figure 10—Cold air drainage effect for the Olympic Peninsula. The map displays the change 

in temperature owing to cold air drainage, expressed as the relative difference in meters of 

elevation.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

21

as elevation plus the latitude and longitude (thus the X, Y, Z 

coordinates) of each geographic 90-m pixel.

The boundary line between the upper boundary of one 

VZ and the lower boundary of another was determined by 

the fit of an equation to achieve a minimum of total errors, 

with a balance of errors above and below the boundary line. 

The boundary equation was typically represented by a linear 

or sine function; however, nonlinear or quadratic functions 

were permitted. In this process, sets of field plots representing 

pairs of adjacent vegetation zones for a small geographical 

area were graphed for elevation versus each environmental 

variable one at a time, starting with PSL. The example in 

figure 4 shows nine PSFZ plots below the lower boundary 

(LBPSFZ) line and 10 WHZ plots above for a total error of 

19 out of 160 total WHZ and PSFZ plots (12 percent error). 

The total number of errors of 19 was the minimum value 

attainable by iteration; thus the errors above and below the 

boundary line were unequal but the minimum attainable. 

At the lower boundary mountain hemlock zone (LBMHZ) 

there were 5 MHZ plots below the upper line (fig. 4) and 3 

PSFZ plots above for a total of 8 out of 72 plots (11 percent 

error). Thus the iteration could not find a solution of eight 

or fewer plots where the error plots above equaled the error 

plots below the boundary line. 

As each variable was added to the model, a new term 

was added resulting in a new value for “A”, “B”, and an 

x+1 number of variables. Each new term was an iteratively 

fit linear or quadratic expression of the relationships of the 

new variable and the elevation of the intercept, or it was a 

complex polynomial reflecting the interaction between it and 

another variable. 

Accuracy assessment of VZ prediction for the ecoplots 

was initially accomplished by comparing the output VZ 

model grid with the VZ call for each plot in the data set. 

Accuracy is calculated by comparing the predicted VZ for 

each plot location with the actual VZ call for each plot. This

same procedure was used for the check plots and ONP plot 

sets when the model was expanded to the Olympic Peninsula. These spatial accuracy assessments are reported as 

percentage accuracy of VZ prediction for the plot data sets.

Model Validation

Two independent sets of plot data were used to validate the 

map output of the VZ model. These FIA and CVS plot data 

were used as a final and unbiased measure of the accuracy 

of the VZ model, as they were installed on the national systematic grid and the data collected were independent of the 

VZ model development. Plots were excluded from the original data sets if they did not have an accurate GPS location 

or could not be assigned a VZ because of stand condition, 

disturbance, or stand age. Of the 132 FIA plots available, 71 

plots off national forest land were usable in this test, and 84 

of the 140 CVS plots on the Olympic National Forest could 

be used for the validation of the PNV model. These 84 CVS 

plots plus 71 FIA plots were compared to the VZ model 

grid for the final spatial accuracy assessment of the model. 

Vegetation zone calls for the plots were compared with the 

mapped VZ for each plot location to determine the percentage accuracy of the VZ prediction.

Results

New Environmental Variables

Four new environmental variables were developed for this 

study: (1) precipitation at sea level (fig. 6), (2) fog effect (fig. 

7), (3) mean annual temperature at sea level (fig. 8), and (4) 

cold air drainage effect (fig. 10). Precipitation at sea level is 

a calculation of the spatial pattern of precipitation with the 

effects of elevation removed. In a practical sense, it depicts 

the precipitation regime for the landscape. The highest PSL 

occurs on the windward or southwestern side of the Olympic Mountains, and becomes successively less going from 

southwest to northeast across the peninsula with an apparent 

cumulative rain-shadow effect behind each ridge. Precipitation at sea level is also less along the coast than farther 

inland toward the mountains. Similarly, MATSL represents 

the spatial pattern of MAT with the effects of elevation 

removed. It depicts the temperature regime for an area. 

The pattern of MATSL on the peninsula shows the lowest 

values at the northwestern corner of the peninsula and along 

the land/saltwater interface. It becomes warmer at greater 

distances from the saltwater, and reaches the highest value 

in the southeastern part of the peninsula farthest from major 
GENERAL TECHNICAL REPORT PNW-GTR-841

22

bodies of saltwater (i.e., Pacific Ocean and Strait of Juan de 

Fuca). Fog effect represents an interpretation of the combined effects of the coastal maritime environment on the 

growth and distribution of plants. These effects are understood to be the results of fog, tree drip, humidity, and cloudiness affecting the water relations of plants. They simulate 

the addition of water in the ecosystem. Cold air drainage 

effect was developed to help explain the downward extent 

of species and vegetation in corridors and basins where cold 

air is believed to flow and accumulate. Its effects are local, 

but when present, CAD effect can have a profound effect on 

the extension of the lower boundaries of most VZs to lower 

elevation. These new environmental variables are presented 

as hypotheses regarding their predictive ability, derivation, 

and conceptual validity. They are not directly measurable, 

nor are they direct derivatives from weather station records; 

therefore, they represent a degree of interpretation. 

Model Building

In the stepwise addition of variables, the base model was 

established and then a new variable was added, one at a 

time, until a final model was achieved. The output of the 

process is summarized in table 1. The base model was the 

simplest version of the model possible, where the entire 

study area was modeled as the most common VZ identified 

from the field data. This was the WHZ, which was represented by 53.1 percent of the USFS ecoplots on the Olympic 

National Forest. Thus if the entire area was modeled as the 

WHZ, it would yield a model and a map with 53.1 percent 

accuracy for the ecoplots. Such a case might represent the 

resolution of a global scale map of PNV with each pixel representing an area of 2 by 3 degrees of latitude and longitude. 

Any model not achieving this level of accuracy was automatically rejected. 

The next level of the model considered only the variable 

elevation, where Yelev corresponded to an intercept value 

for each VZ boundary in a lookup table (thus Yelev = Xelev). 

The elevation values for the boundary equation intercepts 

were lower boundary alpine zone (LBAlpZ) = 1835 m; lower 

boundary subalpine parkland zone (LBPklZ) = 1447 m; 

LBMHZ = 1148 m; LBPSFZ and upper boundary WHZ = 

Table 1—Accuracy assessment for model stages, variables, and plot sets 

677 m. These intercept values were derived by first averaging 

Correct Correct Correct Correct

Correct Accuracy check ONP CVS FIA plots

Number of Number Variable ecoplots (of added for plots (of plots (of plots (of (of 71

Model stagea variables of VZ addeda 1,497 plots) ecoplots 447 plots) 1,051 plots) 84 plots) plots)

Percent

1 - All as single VZ (WHZ) 0 1 53.1 53.12 - Elevation 1 5 Elevation 57.6 4.5 51.5 52.8

3 - Elev and PSL 2 8 PSL 69.6 12.0 69.7 54.6

4 - Elev, PSL, TM 3 8 TM 71.5 1.9 73.9 58.9

5 - Elev, PSL, TM, FOG 4 8 FOG 74.4 2.9 80.9 58.2

6 - Final model 8 8 CAD, aspect, 76.4 2.0 81.3 69.1 82.1 71.8

SW, MATSL a VZ is vegetation zone, ONP is Olympic National Park, CVS is Current Vegetation Survey, FIA is Forest Inventory and Analysis, WHZ is western hemlock zone, Elev is elevation, PSL is precipitation at sea 

level, TM is topographic moisture, FOG is fog effect, CAD is cold air drainage effect, SW is potential shortwave solar radiation, and MATSL is mean annual temperature at sea level.Notes: The first model stage (single VZ model) illustrates the prediction accuracy of the VZ for the ecoplots if the entire study area were mapped as the WHZ and as a single pixel. In this case, 53.1 percent of the ecoplots would have the VZ correctly identified and mapped. With expansion to the 90- by 90-m pixel resolution and the addition of elevation as a lookup variable, the accuracy of ecoplot prediction rises to 57.6 percent. By adding PSL to elevation and using the simple model Y = intercept (elev) + slope x PSL, the accuracy of prediction of the ecoplots rises to 69.6 percent. In the final model stage with seven variables plus the intercept (as elevation), the prediction accuracy for ecoplots rises to 76.4 percent and is 82.1 percent for an independent set of CVS plots on the Olympic National Forest. The average accuracy for the 

combined CVS and FIA validation data sets is 77.4 percent.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

23

the elevation for each VZ (table 2), then averaging the elevation of each adjacent pair to get the initial lower or upper 

boundary for each VZ, and then using the minimum total of 

errors and balance of errors above and below the boundary 

curve, as described previously. Using this approach to define 

the boundary between each of the five dominant VZ, the 

model accuracy for the ecoplots was 57.6 percent (table 1). 

The Sitka spruce, Douglas-fir, and subalpine fir zones were 

not included in this model because they did not have discrete 

boundary elevations; these VZs were limited to areas of 

high fog effect and PSL (SSZ) or to dry TMs or low PSL 

(DFZ and SAFZ). Thus only five of the eight VZs had any 

area in this model that used elevation as the sole variable. 

Subsequent levels of the model were built with a 

stepwise addition of variables. The next level of the model 

added the variable PSL. This yielded a model for eight VZs 

with 69.6 percent accuracy when the model was compared 

to the ecoplots used to build the model, an improvement of 

12 percent over the elevation-only model (table 1). Furthermore, this model included all eight VZs. This large jump 

in accuracy indicated that besides the effect of elevation in 

controlling the intercept values of the boundary equation 

model, PSL added the greatest contribution to the model 

predictability. The next model with PSL and TM as independent variables yielded an accuracy of 71.5 percent. The next 

model added fog effect. This model with only PSL, TM, and 

fog effect as independent variables had an accuracy of 74.4 

percent for the ecoplots, only 2 percent lower than the full 

model, indicating that most of the ecological information 

beyond elevation was contained in just these three independent variables. Finally, CAD effect, aspect, MATSL, and 

shortwave solar radiation were added to the model. In this 

final model, these four latter variables added a combined 

increase of 2 percent in model accuracy for a total accuracy 

of 76.4 percent for the ecoplots used to build the model 

(table 1). 

The single most important independent variable in this 

model was PSL. This simple model Yelev = A + B × (PSL), 

correctly predicted the VZ for 69.6 percent of the ecoplots 

(table 1). The importance of this variable can be seen in the 

relationships shown between PSL, aspect, and elevation and 

the lower boundary of the PSFZ (fig. 11). The lower boundary of the PSFZ changes about 200 m in elevation for each 

class of PSL. Across the widest range of environments from 

the driest area of the peninsula to the wettest, this amounts 

to a total of about 1000 m of elevation difference for the 

lower boundary of the PSFZ, all other factors held constant.

Boundary Equations

The VZ model presented here used boundary equations 

to determine the elevation for each VZ boundary for each 

pixel in the study area. The general form of the final model 

algorithm is given as equation 3, where elevation of the VZ 

boundary is the dependent variable. Table 3 presents the 

specific equations and coefficients for each VZ boundary on 

the Olympic Peninsula. As an example, the following (eqn. 

Table 2—Area and mean elevation of vegetation zones on the 

Olympic Peninsula, Washington

Vegetation 90-m Mean Percentage

zone pixels Area elevation of area 

 Count Hectares Meters Percent

Alpine 15,255 12,357 1,923.2 0.8

Subalpine parkland 67,304 54,518 1,513.6 3.7

Subalpine fir 17,557 14,222 1,542.9 1.0

Mountain hemlock 88,740 71,882 1,230.0 4.8

Pacific silver fir 305,348 247,341 738.3 16.6

Douglas-fir 9,247 7,490 743.6 0.5

Western hemlock 929,336 752,789 265.5 50.5

Sitka spruce 407,608 330,174 78.9 22.1

Sum/mean 1,840,395 1,490,773 411.3 100.0
GENERAL TECHNICAL REPORT PNW-GTR-841

24

4) is the specific model algorithm for the lower boundary of 

the PSFZ, from equation 3 and table 3 E5:

 LBPSFZelev = A + B × (PSL + FOG) + C × TM + D 

 x MATSL + G × FOG + H × CAD + {P x [E + F × 

 (PSL + FOG)] × SIN (ASPECT - 120) + [(1 - P) × 

 E × (SW - 14930) / 3380.4]} (4)

Where

A = 1750; B = -0.300; C = -76.2; D = 10.97; E = 97.5; F = 

-0.0132; G = -0.180; H = -0.800, and P = 1.0; and where LBPSFZelev is the elevation in meters of the lower boundary of 

the PSFZ for any particular 90-m pixel in the study area; the 

corresponding values for each environmental variable for 

the pixel: where (PSL + FOG) is precipitation at sea level in 

millimeters of precipitation plus fog effect in millimeters of 

water equivalent for the contribution of fog effect; TM is the 

topographic moisture for the same pixel; MATSL is the mean 

annual temperature at sea level in degrees C for the same 

pixel; CAD is the cold air drainage effect in meters of elevation change; and SIN (ASPECT-120) is the transformation 

of aspect to set the maximum aspect effect at 210 degrees 

and at 30 degrees; and SW is the shortwave solar radiation in 

watts per square meter, for each 90-m pixel in the landscape calculated using the model application of Kumar et al. 

(1997). The values 14930 and 3380.4 are used to calibrate 

the solar radiation values to relative values of the same scale 

as the sine of the aspect.

Map of Potential Vegetation Zones for the Olympic 

Peninsula

The principal output of the model was a spatial grid (map) of 

the eight VZs of the Olympic Peninsula (fig. 12). The model 

algorithm used to generate this grid computed the boundary elevations for each VZ in the study area by executing 

each of the equations (e.g., table 3) for each of the 1,840,395 

90-m pixels in the study area. It then compared the calculated minimum and maximum elevation of each VZ with the 

elevation of the pixel and selected the VZ whose upper and 

lower boundary bracketed the elevation of the pixel. The 

result was a VZ grid (map) with a 76.4 percent accuracy 

Figure 11—Elevation of the lower boundary of the Pacific silver fir zone by aspect, for seven classes 

of precipitation at sea level (PSL). The five middle classes each represent 508 mm of PSL, and the 

highest and lowest classes are represented by the actual mean PSL for the class.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

25

compared to the USFS ecoplots used to build the model 

(table 1). The predicted VZ accuracy of the check plots used 

to expand plot coverage off the national forest was 81.3 

percent. Accuracy of the combined plot sets from Olympic 

National Park was 69.1 percent (table 1). The final assessment was a comparison with independent plot sets from 

the USFS FIA/CVS program. These 84 CVS plots (national 

forest land) were correctly predicted 82.1 percent of the 

Table 3—Final model coefficients and equations for vegetation zone (VZ) boundary equations for the Olympic 

Peninsula, Washington

 Model coefficients

VZ boundary name Eqna A B C D E F G H P

Lower boundary

alpine zone E1 1834 -0.024 -1.52 1.65 6.10 0.3

Lower boundary

subalpine parkland zone E1 2847 -0.540 -6.10 10.97 18.3 0.3

Upper boundary

subalpine fir zone E1 1911 -0.372 -16.8 5.49 18.3 0.5

Lower boundary

subalpine fir zone E1 1402 0.168 -13.7 0.0 -42.7 0.5

Upper boundary

mountain hemlock zone E4 1908 -9.144 -616257 541.93 -274 2.440 0.5

Lower boundary

mountain hemlock zone E2 1677 -0.216 -22.9 8.78 19.8 -0.060 -0.250 0.3

Lower boundary

Pacific silver fir zone E5 1750 -0.300 -76.2 10.97 97.5 -0.0132 -0.180 -0.800 1.0

Upper boundary

Douglas-fir zone E3 1954 -1.032 -91.4 173.74 0.3

Lower boundary

Douglas-fir zone E2 -1440 0.960 -106.7 -17.56 -183 0.240 2.100 0.3

Upper boundary

Sitka spruce zone E2 -647 0.204 45.7 -29.63 -10.7 2.640 1.000 0.5 a Equation number. Equations follow the form as given in equation 3 in text, where the dependent variable (the left side of the equation) is in 

units of elevation. Each different equation represents the lower and/or upper boundary for one or more vegetation zones.

E1: Yelev = A + B × (PSL + FOG) + C × TM + D × MATSL + {P × E × SIN(ASPECT - 120) + [(1 - P) × E × (SW - 14930) / 3380.4]}

E2: Yelev = A + B × (PSL + FOG) + C × TM + D × MATSL + G × FOG + H × CAD + {P × E × SIN(ASPECT - 120) + [(1 - P) x E × (SW - 14930) / 

3380.4]}

E3: Yelev = A + B × (PSL + FOG) + C × TM + {P × D × SIN(ASPECT - 120) + [(1 - P) × D × (SW - 14930) / 3380.4]}

E4: Yelev = A + B × TM + (C + D × TM 3

) / (PSL + FOG) + {(P × E + G × TM 2) × SIN(ASPECT - 120) + [(1 - P) x E × (SW - 14930) / 3380.4]}

E5: Yelev = A + B × (PSL + FOG) + C × TM + D × MATSL + G × FOG + H × CAD + {P × [E + F × (PSL + FOG)] × SIN(ASPECT - 120) + [(1 - P) 

× E × (SW - 14930) / 3380.4]}

where capital letters refer to coefficients given in above table, and variables identified as: 

PSL + FOG = precipitation at sea level plus fog effect (mm); TM = topographic moisture; MATSL = mean annual temperature at sea level (°C); 

FOG = fog effect (mm); CAD = cold air drainage effect (m); SW = potential shortwave solar radiation (watts per square meter); coefficient “P” is 

proportion of “aspect” effect (“E”) attributed to aspect and “1 - P” is proportion of “aspect” effect (“E”) attributed to potential shortwave solar 

radiation. 

time, and the 71 FIA plots (non-national-forest land) were 

correctly predicted 71.8 percent of the time (table 1). 

The most abundant VZ was the western hemlock zone 

covering 752 789 ha, or about 50.5 percent, well-distributed 

throughout the study area. The Sitka spruce zone was 

second most common with 330 174 ha or about 22 percent of 

the area. It is abundant in the far western Olympic Peninsula. The least common VZs were the alpine zone and the 

Douglas-fir zone, each with less than 1 percent of the study 
GENERAL TECHNICAL REPORT PNW-GTR-841

26

area (table 2). The alpine zone is mostly snow and rock in 

the central and western peninsula, but is represented by 

often well-developed alpine vegetation in the northeastern 

peninsula. The DFZ is mostly restricted to the dry northeastern peninsula, and is usually restricted to south aspects 

on low topographic moisture sites. The PSFZ occurs mainly 

in the central peninsula, at mid elevations and in areas with 

persistent winter snowpack. It is the third most common VZ 

on the peninsula. The MHZ occurs above the PSFZ for most

of the peninsula, and below the subalpine parkland zone. 

The SAFZ is mostly restricted to the dry northeast part of 

the peninsula, replacing the PSFZ or the MHZ there. The 

subalpine parkland zone occupies the tops of most higher 

elevation ridges but occurs below the alpine zone. It is usually a mosaic of tree islands and shrub or herb meadows, 

marked at its lower limit by closed or nearly closed forest 

and at its upper boundary by krummholz, where tree species 

occur as shrubs.

Figure 12—Potential natural vegetation zones of the Olympic Peninsula. The map is the output 

of the model algorithms for the eight potential vegetation zones on the peninsula.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

27

Discussion 

The development of a landscape vegetation model based on 

boundary equations represents a unique approach, as far as 

the authors know. Four of the model variables are original 

(PSL [fig. 6], fog effect [fig. 7], MATSL [fig. 8] and CAD 

effect [fig. 10]). In addition, this calculation of temperature 

lapse rate (fig. 9) and the scaling of topographic moisture 

were also new. The VZ map of the Olympic Peninsula 

(fig. 12) is the finest scale and most accurate map of this 

type, to date, for the Pacific Northwest. 

The original purpose of this project was to develop a 

reliable and accurate map of the PNV of the Olympic 

National Forest through a combination of field reconnaissance and landscape modeling. It began in 1982 drawing on 

experience from two earlier landscape-modeling projects 

in Olympic National Park (Cibula and Nyquist 1987) and in 

Utah (Davis 1980). At that time, computer technology did 

not allow the full expression of the envisioned model.

Without sufficient computer resources at the time, but 

with a model in mind, the alternative was to test the model 

using hand computations and apply the results to topographic maps. Similarity or nearest-neighbor approaches did 

not lend themselves to this field-based approach. However, 

a boundary equation approach did, although nothing was 

found in the literature or from knowledgeable sources to 

support it. Existing variables were scant and consisted of 

slope, aspect, elevation, and an early precipitation map. 

Weather stations on the Olympic Peninsula were few and 

were located mostly at low elevations near saltwater. None 

of these variables provided sufficient information to do a 

reliable map of vegetation zones for this area. Thus a 

new approach was initiated to test a method to distinguish 

between different types of vegetation, and to develop new 

variables that might have better predictability and better 

represent the operational environment of these plant communities. A similar search for variables more explanatory of 

the operational environment was also advocated by Lookingbill and Urban (2005).

We pursued a better identification of the variables that 

compose the operational environment. At the same time, we 

also strove for a better quantification of the mathematical 

relationships between highly correlative environmental 

factors and plant distribution and abundance. That we came 

up with the precipitation regime being most correlated with 

vegetation distribution patterns should not be surprising. 

What is more informative, however, is how we were able 

to build a predictive model with a broad landscape variable 

(PSL), filter that factor out and then analyze the “errors.” 

The subsequent patterns that emerged revealed more about 

the factors driving the vegetation distribution. This led us 

to evaluate the residual pattern at each step, and then begin 

a search for a variable or effect that could be quantified. We 

believe this process provided better insights to the operational environment than we would have been able to derive 

by starting with a set of “existing” variables and performing 

a statistical analysis, weighting each variable similarly, and 

looking for a statistical pattern of correlation.

This model used boundary algorithms to predict the 

environmental limits of the eight VZs modeled in the study 

area. The use of “boundary” in this context refers to an 

ecological node that is used to distinguish between the eight 

VZs in this study. It does not refer to vegetation boundaries 

in the sense of physical ecotones (Carter et al. 1994, Choesin 

and Boerner 2002), ecological edges (Fortin 1994), vegetation boundaries (Glavac et al. 1992), or discontinuities along 

ecological edges (Ludwig and Cornelius 1987). Rather it 

refers explicitly to algorithm-defined boundaries between 

units of PNV. Such boundary equations are used to identify points on the ground and in ecological hyperspace that 

represent these nodes. In addition, the development and use 

of the boundary equations (table 3) provides a new opportunity to examine the relationships between environmental 

variables across geographic and ecological space, the interaction of environmental variables and vegetation, as well as 

changes over time relative to climate change. 

Associations between the different VZ boundaries and 

individual environmental variables and their interactions 

revealed patterns that help to understand the ecological links 

between vegetation and the landscape. One aspect of vegetation and its relationship to the environment that is often 

overlooked in similar studies is that the relative importance 

of different variables changes with different units (i.e., VZs) 
GENERAL TECHNICAL REPORT PNW-GTR-841

28

of the vegetation continuum and in context with other variables (table 3). Observations from these data showed PSL 

was strongly correlated with the distribution of these eight 

VZs except for the alpine zone, which was most strongly 

correlated with elevation (i.e., temperature). The SSZ was 

most correlated with the fog effect, and somewhat with 

MATSL and aspect. The DFZ was defined by both an upper 

and lower boundary equation and was also strongly affected 

by aspect, TM, and MATSL at both boundaries, and by a 

negative CAD effect (i.e., the DFZ avoided areas with high 

CAD effect). The PSFZ was strongly affected by TM, CAD, 

and fog effect as well as PSL. It was also one of only two 

VZs that showed a strong aspect effect, plus it showed a 

strong interactive term between aspect and PSL.

Modeling the boundaries of VZs where the zone is 

continuous around the 360-degree aspect cone is done with 

continuous variables where the upper boundary of one VZ 

is the lower boundary of another. However, when the VZ is 

at the tail of its ecological distribution, it may only occur 

on the drier, southwesterly aspects if at the northerly or wet 

end of its range; conversely, it may only occur on the coolest 

northerly aspect or highest TM at its southerly, warm or 

dry end of its range. In these cases, the boundaries of the 

VZ are described by both a lower and an upper VZ equation. 

Mathematically the wedge formed by the lower and upper 

boundary curve intrudes into the space otherwise occupied 

by another more geographically zonal VZ. This was true for 

the Sitka spruce zone, Douglas-fir zone, and the subalpine 

fir zone in the study area.

Precipitation at sea level was the most predictive variable in this model, indicating the importance of this index 

of the precipitation regime in determining the distribution 

patterns of the VZs in the Olympic Peninsula. The lower 

boundary of the PSFZ (upper boundary of the WHZ) varies 

from the southwest Olympics, where it is lowest, to northeast, where it is absent on some aspects and TMs. Across 

this geographical and topographical gradient, the elevation of the lower boundary increases with decreasing PSL 

(fig. 11) and fog effect and increasing MATSL, and shifts 

upward approximately 100 m on southwesterly aspects. The 

breadth of the PSFZ changes from wetter to drier areas. 

This is shown when comparing figure 5B (an area of high 

PSL) to figure 5C (an area of relatively low PSL). The upper 

boundary of the PSFZ shifts across this PSL gradient, but 

at a lower rate than the lower boundary, resulting in a contraction of the width of the zone from wetter to drier areas, 

approaching the extinction point at about 742 mm PSL, 

where the PSFZ is restricted to high TM.

This model confirms the relationship most often found 

in other studies, that some measure of moisture availability 

or potential moisture loss (Barton 1994, Bridge and Johnson 

2000, Gagnon and Bradfield 1987, Rehfeldt et al. 2006) is 

the variable that is most strongly associated with the distribution of species and communities across the landscape. 

This should cause modelers of global temperature change to 

reconsider the relative effects of temperature versus precipitation in causing major distribution shifts or extinction of 

plant species in the future.

One might expect that the upward shift of VZ boundaries would be proportional to the temperature lapse rate for 

an area, but that was not found. The actual rate of shift upward in elevation is usually less than the temperature lapse 

rate for the area. The explanation appears to be the strong 

correlation between the environmental variables and geography, and the complex interactions between temperature and 

moisture with elevation, each compensating for the other in 

as yet unknown terms. 

The two new variables of PSL and MATSL can be 

regarded as descriptors of the precipitation regime and 

the temperature regime, respectively. As such, they are 

broad-scale variables that interact with the finer scale site 

variables, such as aspect, TM, and CAD effect. As broadscale variables, they are best associated with broad patterns 

of vegetation. These are patterns that would be seen from 

a distance or at a small mapping scale. Within this broad 

context there is considerable variation owing to a different 

set of factors that affect micro- or meso-scale variations in 

light, heat, soil moisture, or nutrients. Using the variables 

PSL and MATSL early in the model-building process is 

not only statistically logical but ecologically logical. The 

effects of microsite variables are easier to understand when 

viewed in the context of the broad-scale factors such as 
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

29

PSL and MATSL. Using PSL and MATSL as indices of 

the precipitation regime and the temperature regime also 

helps with logical consistency. Both variables portray the 

broad precipitation or temperature regime with the effects of 

elevation removed. As elevation is the dependent variable in 

this spatial model, it is hard to justify using variables such 

as total annual precipitation or mean annual temperature as 

independent variables, because both are strongly correlated 

with elevation. Having elevation also strongly represented in 

one or more of the independent variables is circular logic. In 

addition, maps of total annual precipitation and mean annual 

temperature are also models with their own set of errors, 

especially in an area with so few weather stations and such 

extreme topography.

One way to visualize the effects of PSL over a landscape is to imagine several different areas with the same 

PSL. As an example, the PSL 2286 mm isoline (fig. 6) represents an area where the pattern of vegetation is similar for 

the same combination of elevation, aspect, TM, and CAD 

effect. But this pattern is distinctly different in drier or wetter PSL areas. Likewise, the MATSL variable represents a 

temperature regime across the landscape, where for a given 

value of MATSL, the pattern of vegetation relative to the site 

variables is relatively consistent. A number of hypothetical 

cone-shaped mountains scattered throughout the area of PSL 

2286 mm and MATSL 10.5 °C would have similar patterns 

of VZs from bottom to top. Likewise a set of mountains in 

the area of PSL 1143 mm and MATSL 10.5 °C would have 

similar species and patterns of vegetation, but the vegetation distribution for the mountain in PSL 2286 mm would be 

distinctly different from those in PSL 1143 mm.

Model validation for the eight VZs of the 1.5-millionha study area using independent plot data was 77.4 percent. 

This was the average of the prediction accuracy of the CVS 

and FIA independent plot sets (table 1). Too few studies in 

the literature use an independent and statistically unbiased 

validation sample to test model accuracy (Reynolds and 

Chung 1986). 

No directly comparable published studies were found. 

However, for a model of current species and structural attributes at a 25-m resolution, Ohmann and Gregory (2002) 

reported a predictive accuracy of 45 percent for a reserved 

set of 25 percent of the original plots. Rehfeldt et al. (2006) 

compared a model assignment of plots to a preexisting polygon map of 25 biotic communities (Brown et al. 1998) in the 

Western United States. They achieved a classification error

of 19 percent when comparing precisely located plots to a 

polygon map of unknown accuracy and geographic registration. Lookingbill and Urban (2005) compared two different 

sets of environmental variables: one a “proxy model” and 

another with more plant-relevant variables, but with a nonrandomly located validation sample (n = 11) to four forest 

community types in western Oregon. They found 73 percent 

accuracy for the “plant-relevant” variables but only 18 percent for the proxy model variables.

This model is easily verified in the field by having map 

in hand. In addition, by knowing the PSL, MATSL, aspect, 

and TM, one can predict with relative confidence the VZ for 

any particular site on the landscape. For the roughly one out 

of five times the model will not correctly identify the fieldidentified VZ (table 2), one can usually deduce one of four 

possible sources of the error. First, the site is very near the 

boundary of the VZ in question, and the coarseness of the 

independent variables or the influence of other micro-site 

variables may explain the difference; two, the community 

might be in an early successional stage or disturbed condition, and it might not be possible to identify the VZ in the 

field; three, the site is affected by a variable not recognized 

in the model, possibly owing to lack of GIS data, such as a 

spring, seep, wind-affected redistribution of precipitation, 

or an unusual or unmapped soil or geological element; or 

fourth, the location of the site is not known with the same 

precision as the GIS-based map. Many of the plots used in 

developing the model were installed before the widespread 

use of GPS and therefore their locations are known with less 

precision than later plots. This locational error is believed to 

account for about half of the statistical “error” of the model. 

The qualitative assignment of VZ to plots, especially those 

done “after the fact,” may be another source of error. When 

the model does not correctly identify the apparent VZ in the 

field and the location of the site is known with a high degree 

of precision, it is more interesting than when it does, as it 
GENERAL TECHNICAL REPORT PNW-GTR-841

30

encourages the pursuit of new ecological information or 

definition of a new environmental variable.

This modeling approach was developed for the broadscale classification of PNV, i.e., the VZ. However it is believed to be equally applicable to finer scale vegetation units 

(i.e., plant association groups [PAG] and plant association). 

In addition, the concepts behind the model presented here 

could also be used to describe

the patterns of distribution of PNV as a noncategorical 

variable. As such, PNV can also be understood in terms of 

a scale or gradient of potential vegetation conditions that 

are not encumbered by artificial classes. Within each VZ 

described here by ecological-mathematical boundaries, there 

is continuous vegetational compositional change associated 

with more or less continuous environmental variation. This 

gradient approach is as much a part of the concept of potential vegetation as the typal-classification approach. Each 

has different uses and values in understanding the nature of 

vegetation or the application of concepts or stratifications to 

natural resource management.

The described patterns of potential VZs provide a broad 

vegetational stratification for addressing ecological or management questions, as the VZ and the broad classes of PSL 

and MATSL provide a gauge of relative ecological similarity. This stratification can provide context for comparing 

results of successional studies. The concept of succession 

describes the natural and inevitable changes that occur in 

plant communities over time. These changes are strongly 

affected by the starting condition and the environment of 

the site and are driven by competition for space and resources. When studying the effects of succession on community 

biomass, growth, diversity, or function, the factor most 

important in controlling these effects is the environment. 

The environment is strongly reflected in the potential vegetation for a site. Thus, knowing the PNV for any stand or 

community being studied helps in stratifying the study and 

in interpreting the results. 

Maps of potential vegetation and environmental variables can provide the foundation for scientific and ecological 

studies and for addressing resource management questions. 

Potential VZs or finer levels of classification of PNV can 

be used to stratify studies of silvicultural treatment, stand 

growth response to disturbance, carbon sequestration, and 

fire resistance or probability. It is the opinion of the authors 

that too few silviculture studies of thinning, weeding, or 

regeneration inform the users regarding where the results 

of a study are applicable. A thinning study of a Douglas-fir 

stand in the grand fir zone in the southern Willamette Valley can hardly be applicable to a similarly aged Douglas-fir 

stand in the moist western hemlock zone in northern Puget 

Sound; yet, such applications have been made, with unfortunate consequences. Stand growth studies are often stratified 

by site index, yet knowing the VZ, plant association group, 

or plant association can inform the user regarding potential 

stockability differences or regeneration potential not associated with site index alone. Recent and unpublished stand 

growth/carbon sequestration studies show that some previously unrecognized stand growth patterns can be stratified 

by units of PNV. Fire history studies of western Washington 

show that historical fire patterns are strongly associated 

with current patterns of PNV. The VZ and PSL data have 

been used to stratify the peninsula into three fire occurrence 

zones: the wet, western part having little evidence of historical fire; the dry northeast part characterized by relatively 

young forests with recurrent forest fires over the last 300 

years; and the middle area, where many forests date from 

either about 1300, 1500, or 1700 a.d. The model of potential 

vegetation presented here represents the current climate. 

However, by varying PSL or MATSL one can use the PNV 

model to predict how PNV could change in the future owing 

to climate change or how it might have been different in 

the past. The PNV model has also been used to predict and 

validate the potential habitat and range for over 10 plant and 

animal species, including species considered rare and at risk 

(Lesher 2005).

This study was a quest for a new conceptual model for 

understanding how vegetation is distributed across a landscape. This model included developing and validating the 

concept of boundary equations, development of several new 

environmental variables, and production of a high-resolution 

grid of potential VZs. This grid map and these variables 

permit new perspectives on the patterns of vegetation and 
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

31

their links to the environment, plus possible uses for studying the consequences of climate change, plant community 

monitoring, and ecological analysis. This project has been 

extended across the rest of Washington and Oregon, where 

the breadth of major environmental variables provides the 

opportunity for exploring nonlinear relationships and factor 

interactions across a regional landscape.

Acknowledgments

The authors wish to acknowledge the assistance of many who 

have helped in this project, including Michael Fitzgibbon of 

ESRI, who helped in converting the original mathematical 

and ecological relationships to Arc macro language. Michael 

was very helpful in making this model application take form 

and function. Others who have helped in writing code for 

the PNV model application include Greg Dillon, Aimee 

Lesieutre, Bill Wettengel, Andy Wilson and Steve Zubalik. The USFS and Bureau of Land Management (BLM) 

personnel who have helped in their support and funding 

include Terry Brumley, Tom Demeo, Nancy Molina, Fred 

Hall, Miles Hemstrom, Randy Molina, John Teply, Bruce 

Rittenhouse, and Lenny Volland. The authors also wishes to 

acknowledge the significant assistance provided by the reviewers Tom Demeo, Jerry Franklin, Nancy Fredricks, Don 

McKenzie, David L. Peterson, Ed Schreiner, and Catharine 

Thompson. Funding was provided by the USFS and BLM 

through the USFS Pacific Northwest Region Ecology Program and the interagency Survey and Manage Program.

Acronyms

ATU Accumulated temperature units

CAD Cold air drainage

CPC Climax plant community

CVS Current Vegetation Survey (USFS)

DEM Digital elevation model

DFZ Douglas-fir zone

ELEV Elevation

FIA Forest Inventory and Analysis (USFS)

GIS Geographic information system

GPS Global positioning system

LBMHZ Lower boundary mountain hemlock zone

LBPklZ Lower boundary subalpine parkland zone

LBPSFZ Lower boundary Pacific silver fir zone

MAT Mean annual temperature 

MATSL Mean annual temperature at sea level

MHZ Mountain hemlock zone

ONP Olympic National Park

PNV Potential natural vegetation

PSFZ Pacific silver fir zone

PSL Precipitation at sea level

SAFZ Subalpine fir zone

SIN Sine function

SSZ Sitka spruce zone

SW Potential shortwave solar radiation

TAP Total annual precipitation

TM Topographic moisture

VZ (Potential) vegetation zone

WHZ Western hemlock zone

English Equivalents

When you know: Multiply by: To find:

Millimeters (mm) 0.0394 Inches

Centimeters (cm) 0.394 Inches

Meters (m) 3.28 Feet

Kilometers (km) 0.621 Miles

Hectares (ha) 2.47 Acres

Square meters (m2) 10.76 Square feet

Square kilometers (km2) 0.386 Square 

 miles

Degrees Celsius (°C) 1.8C + 32 Degrees 

 Fahrenheit

Watts per square meter 0.0929 Watts per 

 square foot

The vegetation zone model and environmental variables presented here are available as maps in grid format at 

a 90-m pixel resolution. They can be downloaded from the 

Ecoshare Web site—http://ecoshare.info/. 

References

Azevedo, J.; Morgan, D.L. 1974. Fog precipitation in 

coastal California forests. Ecology. 55: 1135–1141.
GENERAL TECHNICAL REPORT PNW-GTR-841

32

Barton, A.M. 1994. Gradient analysis of relationships 

among fire, environment, and vegetation in a southwestern 

USA mountain range. Bulletin of the Torrey Botanical 

Club. 121: 251–265.

Beven, K.J.; Kirky, M.J. 1979. A physically based, variable 

contributing area model of basin hydrology. Hydrological 

Science Bulletin. 24: 43–69.

Bridge, S.R.J.; Johnson, E.A. 2000. Geomorphic principles 

of terrain organization and vegetation gradients. Journal 

of Vegetation Science. 11: 57–70.

Brockway, D.G.; Topik, C.; Hemstrom, M.A.; 

Emmingham, W.H. 1983. Plant association and 

management guide for the Pacific silver fir zone, Gifford 

Pinchot National Forest. Tech. Pap. R6-Ecol-130a-1983. 

Portland, OR: U.S. Department of Agriculture, Forest 

Service, Pacific Northwest Region. 122 p.

Brown, D.E.; Reichenbacher, F.; Franson, S.E. 1998. A 

classification of North American biotic communities. Salt 

Lake City: University of Utah Press. 141 p.

Cajander, A.K. 1926. The theory of forest types. Acta 

Forestalia Fennica. 29: 1–108.

Campbell, S.; Waddell, K.; Gray, A., eds. 2010. 

Washington’s forest resources, 2002–2006. Gen. Tech. 

Rep. PNW-GTR-800. Portland, OR: U.S. Department of 

Agriculture, Forest Service, Pacific Northwest Research 

Station. 189 p.

Carter, V.; Gammon, P.T.; Garrett, M.K. 1994. Ecotone 

dynamics and boundary determination in the Great Dismal 

Swamp. Ecological Applications. 4: 189–203.

Choesin, D.N.; Boerner, R.E.J. 2002. Vegetation boundary 

detection: a comparison of two approaches applied to field 

data. Plant Ecology. 158: 85–96.

Chung, U.; Seo, H.H.; Hwang, K.H.; Hwang, B.S.; Choi, 

J.; Lee, J.T.; Yun, J.I. 2006. Minimum temperature 

mapping over complex terrain by estimating cold 

air accumulation potential. Agricultural and Forest 

Meteorology. 137: 15–24.

Cibula, W.G.; Nyquist, M.O. 1987. Use of topographic 

and climatological models in a geographical data base to 

improve landsat MSS classification for Olympic National 

Park. Photogrammetric Engineering and Remote Sensing. 

53: 67–75.

Daly, C.; Gibson, W.P.; Taylor, G.H.; Johnson, G.L.; 

Pasteris, P. 2002. A knowledge-based approach to the 

statistical mapping of climate. Climate Research. 

22: 99–113.

Daubenmire, R. 1952. Forest vegetation of northern Idaho 

and adjacent Washington, and its bearing on concepts of 

vegetation classification. Ecological Monographs. 22: 

301–330.

Daubenmire, R. 1976. The use of vegetation in assessing 

the productivity of forest lands. The Botanical Review. 

42: 115–143.

Daubenmire, R.; Daubenmire, J.B. 1968. Forest vegetation 

of eastern Washington and northern Idaho. Technical 

Bulletin 60. Pullman, WA: Washington State University, 

Washington Agricultural Experiment Station. 104 p.

Davis, L.S. 1980. Strategy for building a locationspecific multi-purpose information system for wildland 

management. Journal of Forestry. 78: 402–408. 

Davis, R.B. 1966. Spruce-fir forests of the coast of Maine. 

Ecological Monographs. 36: 80–94.

Dawson, T.E. 1998. Fog in the California redwood forest: 

ecosystem inputs and use by plants. Oecologia. 

117: 476–485.

Diaz, N.M.; High, C.T.; Mellen, T.K.; Smith, D.E.; Topik, 

C. 1997. Plant association and management guide for 

the mountain hemlock zone. Tech. Pap. R6-MTH-GPTP-08-95. Portland, OR: U.S. Department of Agriculture, 

Forest Service, Pacific Northwest Region. 101 p. 

Dodwell, A.; Rixon, T.F. 1902. Forest conditions in the 

Olympic Forest Reserve, Washington. Professional 

Paper No. 7, Series H, Forestry 4. Washington, DC: U.S. 

Department of the Interior, Geological Survey. 298 p.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

33

Dymond, C.C.; Johnson, E.A. 2002. Mapping vegetation 

spatial patterns from modeled water, temperature and 

solar radiation gradients. Journal of Photogrammetry and 

Remote Sensing. 57: 69–85.

Environmental Systems Research Institute, Inc. [ESRI]. 

1982–2006. ARC 9.2. Redlands, CA. 

Eyre, F.H. 1980. Forest cover types of the United States and 

Canada. Bethesda, MD: Society of American Foresters. 

148 p. 

Forest and Rangeland Renewable Resources Planning 

Act of 1974 [RPA]; 16 U.S.C. 1601 (note).

Fortin, M-J. 1994. Edge detection algorithms for twodimensional ecological data. Ecology. 75: 956–965.

Franklin, J. 1995. Predictive vegetation mapping: 

geographic modeling of biospatial patterns in relation to 

environmental gradients. Progress in Physical Geography. 

19: 474–499.

Franklin, J.F.; Dyrness, C.T. 1973. Natural vegetation 

of Oregon and Washington. Gen. Tech. Rep. PNW-8. 

Portland, OR: U.S. Department of Agriculture, Forest 

Service, Pacific Northwest Forest and Range Experiment 

Station. 417 p.

Franklin, J.F.; Moir, W.H.; Hemstrom, M.A.; Greene, 

S.E.; Smith, B.G. 1988. The forest communities of 

Mount Rainier National Park. Scientific Monograph 

Series No. 19. Washington, DC: U.S. Department of the 

Interior, National Park Service. 194 p. 

Gagnon, D.; Bradfield, G.E. 1987. Gradient analysis of 

west central Vancouver Island forests. Canadian Journal 

of Botany. 65: 822–833.

Glavac, V.; Grillenberger, C.; Hakes, W.; Ziezold, H. 

1992. On the nature of vegetation boundaries, undisturbed 

flood plain forest communities as an example—a 

contribution to the continuum/discontinuum controversy. 

Vegetatio. 101: 123–144.

Gleason, H.A. 1926. The individualistic concept of the plant 

association. Bulletin of the Torrey Botanical Club. 

52: 7–26.

Harr, R.D. 1982. Fog drip in the Bull Run Municipal 

Watershed, Oregon. Water Resources Bulletin. 

18: 785–789.

Henderson, J.A.; Lesher, R.D.; Peter, D.H.; Shaw, D.C. 

1992. Field guide to the forested plant associations of 

the Mt. Baker-Snoqualmie National Forest. Tech. Pap. 

R6-ECOL-028-91. Portland, OR: U.S. Department of 

Agriculture, Forest Service, Pacific Northwest Region. 

196 p.

Henderson, J.A.; Peter, D.H.; Lesher, R.D.; Shaw, D.C. 

1989. Forested plant associations of the Olympic National 

Forest. Tech. Pap. R6-ECOL-001-88. Portland, OR: 

U.S. Department of Agriculture, Forest Service, Pacific 

Northwest Region. 502 p. http://ecoshare.info/olympic-nfplant-associations/ (November 2010).

Hironaka, M. 1985. Habitat type, range site and community 

type. In: Proceedings, Symposium on the biology of 

Artemisia and Chrysothamnus. Gen. Tech. Rep. INT-200. 

Fort Collins, CO: U.S. Department of Agriculture, Forest 

Service, Rocky Mountain Research Station: 15–18.

Hitchcock, C.L.; Cronquist, A.R. 1973. Flora of the Pacific 

Northwest. Seattle, WA: University of Washington Press. 

730 p.

Hungerford, R.D.; Nemani, R.R.; Running, S.W.; 

Coughlan, J.C. 1989. MTCLIM: a mountain 

microclimate simulation model. Res. Pap. INT-414. 

Ogden, UT: U.S. Department of Agriculture, Forest 

Service, Intermountain Research Station. 52 p.

Jennings, M.; Loucks, O.; Peet, R.; Faber-Langendoen, 

D.; Glenn-Lewin, D.; Grossman, D.; Damman, A.; 

Barbour, M.; Pfister, R.; Walker, M.; Talbot, S.; 

Walker, J.; Hartshorn, G.; Waggoner, G.; Abrams, 

M.; Hill, A.; Roberts, D.; Tart, D.; Rejmanek, M. 

2003. Guidelines for describing associations and alliances 

of the U.S. national vegetation classification. Version 

3.0. Washington, DC: Ecological Society of America, 

Vegetation Classification Panel. 
GENERAL TECHNICAL REPORT PNW-GTR-841

34

Jenson, S.K.; Domingue, J.O. 1988. Extracting topographic 

structure from digital elevation data for geographic 

information system analysis. Photogrammetric 

Engineering and Remote Sensing. 54: 1593–1600. 

Kessel, S.R. 1979. Gradient modeling: resource and fire 

management. New York: Springer-Verlag. 432 p.

Kuchler, A.W. 1964. Manual to accompany the map: 

potential natural vegetation of the conterminous United 

States. Special Publication No. 36. New York: American 

Geographical Society. 116 p.

Kumar, L.; Skidmore, A.K.; Knowles, E. 1997. Modelling 

topographic variation in solar radiation in a GIS 

environment. International Journal of Geographical 

Information Science. 11: 475–497.

Lesher, R.D. 2005. An environmental gradient model 

predicts the spatial distribution of potential habitat for 

Hypogymnia duplicata in the Cascade Mountains of 

northwestern Washington. Seattle, WA: University of 

Washington. 85 p. Ph.D. dissertation.

Lillybridge, T.R.; Kovalchick, B.L.; Williams, C.K.; 

Smith, B.G. 1995. Field guide for forested plant 

associations of the Wenatchee National Forest. Gen. Tech. 

Rep. PNW-GTR-359. Portland, OR: U.S. Department of 

Agriculture, Forest Service, Pacific Northwest Research 

Station. 336 p.

Lookingbill, T.R.; Urban, D.L. 2005. Gradient analysis, the 

next generation: towards more plant-relevant explanatory 

variables. Canadian Journal of Forest Research. 

35: 1744–1753.

Ludwig, J.A.; Cornelius, J.M. 1987. Locating 

discontinuities along ecological gradients. Ecology. 

68: 448–450.

Mason, H.L.; Langenheim, J.L. 1957. Language analysis 

and the concept “environment.” Ecology. 38: 325–340. 

McKenzie, D.; Peterson, D.W.; Peterson, D.L. 2003.

Modeling conifer species distributions in mountain forests 

of Washington State, U.S.A. Forestry Chronicle. 

79: 253–258. 

Mueller-Dombois, D.; Ellenberg, H. 1974. Aims and 

methods of vegetation ecology. New York: John Wiley 

and Sons. 547 p.

Oberlander, G.T. 1956. Summer for precipitation in San 

Francisco Peninsula. Ecology. 37: 851.

Ohmann, J.L.; Gregory, M.J. 2002. Predictive mapping 

of forest composition and structure with direct gradient 

analysis and nearest neighbor imputation in coastal 

Oregon, U.S.A. Canadian Journal of Forest Research. 

32: 725:741.

Ohmann, J.L.; Gregory, M.J.; Spies, T.A. 2007. Influence 

of environment, disturbance and ownership on forest 

vegetation of coastal Oregon. Ecological Applications. 

17: 18–33.

Oksanen, J.; Minchin, P.R. 2002. Continuum theory 

revisited: What shape are species responses along 

ecological gradients? Ecological Modelling. 

157: 119–129.

Pellant, M.; Shaver, P.; Pyke, D.A.; Herrick, J.E. 2005.

Interpreting indicators of rangeland health, version 4. 

Technical Reference 1734-6. BLM/WO/ST-00/001+1734/

RE05. Denver, CO: U.S. Department of the Interior, 

Bureau of Land Management, National Science and 

Technology Center. 122 p. 

Peterson, D.L.; Schreiner, E.G.; Buckingham, N.M. 1997.

Gradients, vegetation and climate: spatial and temporal 

dynamics in the Olympic Mountains, U.S.A. Global 

Ecology and Biogeography Letters. 6: 7–17.

Pfister, R.D.; Arno, S.F. 1980. Classifying forest habitat 

types based on potential climax vegetation. Forest 

Science. 26: 52–70.

Pfister, R.D.; Kovalchik, B.L.; Arno, S.F.; Presby, R.C. 

1977. Forest habitat types of Montana. Gen. Tech. Rep. 

INT-34. Ogden, UT: U.S. Department of Agriculture, 

Forest Service, Intermountain Forest and Range 

Experiment Station. 174 p.
A Landscape Model for Predicting Potential Natural Vegetation of the Olympic Peninsula USA

35

Rehfeldt, G.E.; Crookston, N.L.; Warwell, M.V.; Evans, 

S. 2006. Empirical analyses of plant-climate relationships 

for the Western United States. International Journal of 

Plant Sciences. 167: 1123–1150.

Rehfeldt, G.E.; Ferguson, D.E.; Crookston, N.L. 2008.

Quantifying the abundance of co-occurring conifers along 

inland Northwest (U.S.A) climate gradients. Ecology. 

89: 2127–2139.

Reynolds, J.M.R.; Chung, J. 1986. Regression 

methodology of estimating prediction error. Canadian 

Journal of Forest Research. 16: 931–938.

Robertson, M.P.; Caithness, N.; Villet, M.H. 2001. A PCAbased modelling technique for predicting environmental 

suitability for organisms from presence records. Diversity 

and Distributions. 7: 15–27.

Smith, E.L.; Johnson, P.S.; Ruyle, G.; Smeins, F.; Loper, 

D.; Whetsell, D.; Child, D.; Sims, P.; Smith, R.; 

Volland, L.; Hemstrom, M.; Bainter, E.; Mendenhall, 

A.; Wadman, K.; Franzen, D.; Suthers, M.; 

Willoughby, J.; Habich, N.; Gaven, T.; Haley, J. 1995. 

New concepts for assessment of rangeland condition. 

Journal of Range Management. 48: 271–282.

Stone, R.G. 1936. Fog in the United States and adjacent 

regions. The Geographical Review. 26: 111–134.

Tabor, R. 1975. Guide to the geology of Olympic National 

Park. Seattle, WA: University of Washington Press. 144 p. 

ter Braak, C.J.F.; Prentice, I.C. 1994. A theory of gradient 

analysis. Advances in Ecological Research. 18: 271–313.

Topik, C.; Halverson, N.M.; Brockway, D.G. 1986.

Plant association and management guide for the 

western hemlock zone, Gifford Pinchot National Forest. 

Tech. Pap. R6-ECOL-230A-1986. Portland, OR: U.S. 

Department of Agriculture, Forest Service, Pacific 

Northwest Region. 132 p. 

Tuxen, R. 1956. Die heutige potentielle naturliche 

Vegetation als Gegenstand der Vegetationskartierung. 

Stolzenau-au/Weser. Angewandte Pflanzensoziologie. 

13: 5–42.

U.S. Department of Agriculture, Forest Service. 1998.

Region 6 Inventory and Monitoring System: field 

procedures for the Current Vegetation Survey. Version 

2.03. Unpublished document. On file with: USDA Forest 

Service, Pacific Northwest Region, Forest Inventory and 

Monitoring, P.O. Box 3623, Portland, OR 97208.

Vogelmann, H.W.; Siccama, T.; Leedy, D.; Ovitt, D.C. 

1968. Precipitation from fog moisture in the Green 

Mountains of Vermont. Ecology. 49: 1205–1207.

Whittaker, R.H. 1956. Vegetation of the Great Smoky 

Mountains. Ecological Monographs. 26: 1–80.

Whittaker, R.H. 1960. Vegetation of the Siskiyou 

Mountains, Oregon and California. Ecological 

Monographs. 30: 279–338.

Whittaker, R.H. 1967. Gradient analysis of vegetation. 

Biological Reviews. 49: 207–264. 

Whittaker, R.H.; Niering, W.A. 1965. Vegetation of the 

Santa Catalina Mountains, Arizona: a gradient analysis of 

the south slope. Ecology. 46: 422–452.

Williams, C.K.; Lillybridge, T.R. 1983. Forested plant 

associations of the Okanogan National Forest. Tech. Pap. 

R6-ECOL-132b-1983. Portland, OR: U.S. Department of 

Agriculture, Forest Service, Pacific Northwest Region. 

139 p. 

Winthers, E.; Fallon, D.; Haglund, J.; DeMeo, T.; 

Nowacki, G.; Tart, D.; Ferwerda, M.; Robertson, G.; 

Gallegos, A.; Rorick, A.; Cleland, D.T.; Robbie, W. 

2005. Terrestrial Ecological Unit Inventory technical 

guide. Washington, DC: U.S. Department of Agriculture, 

Forest Service, Washington Office, Ecosystem 

Management Coordination Staff. 245 p.