1.Question asked and data used: 

• What are the spatial patterns of juniper sapling establishment within a treated watershed 14 years after most juniper were removed?
• Two data sets were used for this analysis:
  • Thirty-three belt transects (30 m in length, 3 m in width) were used to assess the density of juniper saplings within the watershed (data in form of excel spreadsheet listing latitude and longitude with number of juniper per transect)
  • An orthomosaic created from UAV imagery of a small area within the watershed (data in form of multispectral raster, with brightness values for red, green, blue, near-infrared, and re-edge wavelengths)
• Watershed map (National Agricultural Imagery Program image from 2016 shown) with belt transects and UAV study area:

2. Approach and tools used: A number of techniques were initially applied to explore the data and to examine which approaches were more (or less) effective.

• Belt transects:
  • Kriging
  • Inverse Distance Weighted (IDW) interpolation
    
  • Spatial autocorrelation (Global Moran’s I)
  • Hot-spot analysis
• Classified raster (UAV orthomosaic):
  • Hot-spot analysis

3. Steps used in analysis:

• Slope and aspect calculation:
  • The slope and aspect tools (Spatial Analyst Tools) were applied to 30 m DEM (from earthexplorer.usgs.gov). This information was noted for hot spots and cold spots but not used for further analysis during this exercise.
  • The extract multi values to points was used to extract the slope, aspect, and elevation value
• Belt transects:
  • Projected data into NAD 1983 UTM Zone 10N
  • The location of each survey was symbolized by color based on the number of juniper found.
  • Kriging (Spatial Analyst Tools and Geostatistical Wizard) was used for initial exploratory analysis to assess general spatial distribution across the watershed
    • The number of juniper per transect used as the z-value field
    • Histogram and QQPlot used assess distribution of the number of juniper per transect for the geostatistical layer
    • For the raster layer, the predicted values were compared to the observed values by extracting the values at each survey location
  • IDW interpolation
    •  Juniper used as input feature
    • Maximum neighbors:15; Minimum neighbors:10
  • Spatial autocorrelation (Global Moran’s I) (Spatial Statistics Tools)
    • Calculated using the number of juniper per transect
    • HTML file created lists Moran’s index, z-score, and p-value
  • Hot Spot Analysis (Getis-Ord Gi) (Spatial Statistics Tools)
    • Conceptualization of spatial relationships: Fixed distance band
    • Distance method: Euclidean distance
    • Default was used for the distance band

• Classified Raster (orthomosaic):
  • Supervised classification (support vector machine) was used to identify juniper within the UAV orthomosaic
  • General procedure for classification: create training samples-> assess using interactive training sample manager and scatterplots->create .ecd file->apply classifier->apply Majority Filter tool
  • A point file was created from pixel clusters identified as juniper within the image
  • Hot Spot analysis (Getis-Ord Gi)
    
    • Conducted to assess areas of concentration of juniper saplings within the sample area
    • Integrate tool used to aggregate juniper data (5 m used for analysis, but other distances were initially tested)
    • Hot Spot analysis tool using aggregated data created in previous step as input layer (other inputs were the same as those used for the belt transect hot spot analysis)

4. Overview of Results:

• Belt Transects
  • Kriging. I used two different approaches to kriging. First, I used the kriging tool under spatial analyst tools and then used the geostatistical wizard to calculate ordinary prediction kriging. The resulting maps using these two approaches were different. In particular, the use of the geostatistical wizard created maps more similar to the IDW while the geostatistical wizard created a map with different contours.
    • Related statistics:
      • minimum:0.05
      • maximum: 6.97
      • mean: 2.83
      • standard deviation: 0.86
    • Spatial Analyst Kriging Map:
    •  
    • Using the export values function, the predicted values of this kriging method at each belt transect site were very similar (within 0.05) to the observed values.

• • • Ordinary Prediction Kriging:
    • 
    • The QQ plot indicates that assumptions of normality may not be met:
    • 
    • The ordinary prediction kriging also tended to overestimate the juniper density based upong the distribution chart:
    • 
      • Related cross-validation/error statistics for the ordinary prediction kriging:
        • Mean -0.0583700520511708
          Root-Mean-Square 2.09329560839956
          Mean Standardized -0.0193126589526469
          Root-Mean-Square Standardized 0.986330022079785
          Average Standard Error 2.09938534113134
  • IDW
    • While general trends were similar to kriging, the size and shape of contours between the methods were different.
    • Mean: 3.08, range: 0.00 to 7.00, standard deviation: 1.16
    • IDW
    • 
  • Spatial autocorrelation (Global Moran’s I)
    • Moran’s I based on juniper per transect: -0.019, with a p-value of 0.92
    • Indicates that the pattern of juniper density is considered random and not significantly different than a random distribution
  • Hot Spot analysis (Getis-Ord Gi)
    • One hot spot found (90% confidence): northwestern aspect (300 degrees) with 8.8% slope
    • One cold spot found (95% confidence): north-northwestern aspect (347 degrees) with 11.5% slope
    • Remaining points were considered insignificant
    • Map of Hot Spot analysis:
    • 

• Classified Raster
  • Hot Spot analysis (Getis-Ord Gi)
    • Within the area covered by the orthomosaic, four hot spots were found:
      • One area with 99% confidence: 11% slope, 325 degree aspect
      • Three areas with 95% confidence: 3% slope, 254 degree aspect; 11% slope, 167 degree aspect; and 8% slope,137 degree aspect
    • UAV Orthomosaic Hot Spot map:
    •  

5. Discussion/Critique of Methods.

Based on the maps produced using the IDW and kriging methods, some general trends in juniper spatial distribution exist within this watershed. For example, we see lower densities in the north-northeastern areas and greater densities of juniper in the southeastern and western areas. However, both the IDW and the kriging raster using the spatial analyst tool produced the characteristic “bulls-eye” pattern often associated with IDW. When compared to the ordinary prediction kriging maps, different interpretations of juniper density in the watershed might be made. Compared to the belt transects data, the kriging created using the spatial analyst tool was similar to the observed values while the ordinary prediction kriging tended to overestimate the distribution. However, more ground data is necessary to determine how accurate these prediction would be in other areas.  One of the biggest takeaways from this is to carefully consider the approach used and the nature of the data (for example, the use of ordinary versus simply kriging, etc.). From a user standpoint, I found the geostatistical wizard the most useful approach — particularly as it made inspecting the statistics and semivariogram very easy. However, in the future I would explore different methods of kriging within this tool.

The spatial autocorrelation and hot spot tools were useful, although the results did not provide much significant information regarding the spatial distribution of juniper in the cases examined here. For future analysis, particularly when examining the relationship between juniper density and other watershed characteristics (e.g., slope, aspect, etc.) these tools may become more important. In the case of the UAV orthomosaic, this provided a “test run” for what will be a larger dataset. The steps taken prior to analysis, particularly creating a point layer from the classified pixel clusters, will be time intensive and may require alternative approaches. At the limited scale of the UAV orthomosaics these hotspots did not provide much useful information, but if extrapolated over a larger area more patterns may be observed.

The distribution of juniper saplings in this case may be associated with a number of factors. Juniper seeds may be spread through birds or wind, resulting in the spread of juniper saplings throughout the watershed. At a much larger scale, if seed sources are less available, this distribution may be more localized. This pattern may also vary with the presence of mature juniper. As previously discussed some patterns emerged in this data (lower densities in the northern areas of the study area, for example) which may be associated with the spatial distribution of other characteristics, such as soil moisture. For example, sources and areas of higher juniper may include areas were we anticipate higher soil moisture such as flatter slopes with northerly aspects. These factors will be assessed further in exercise 2.