<jupyter_start><jupyter_text>Simulating patterned vegetationUsing the model described in the paper by Konings *et al.* (doi :10.1029/2011JG001748), we can define a 50x50 grid, in which each square contains values for surface water, soil water, and plant biomass, and see how the plant biomass patterns evolve as we progress the simulation, given different values of daily rainfall.We can then use the subgraph centrality method described in Mander *et al.* (https://doi.org/10.1098/rsos.160443) to quantify the *connectedness* of the vegetation patterns.<jupyter_code>import os
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline<jupyter_output><empty_output><jupyter_text>In the pyveg package we have a class called PatternGenerator that implements the model:<jupyter_code>from pyveg.src.pattern_generation import PatternGenerator
pg = PatternGenerator()<jupyter_output><empty_output><jupyter_text>Once we have an instance of PatternGenerator, we can print its configuration:<jupyter_code>pg.print_config()<jupyter_output>Current configuration:
======================

m: 50
delta_x: 21
delta_y: 21
diffusion_plant: 0.1
diffusion_soil: 0.1
diffusion_surface: 100
surface_water_frac: 0.1
bare_soil_infiltration: 0.15
grazing_loss: 0.3
soil_water_loss: 0.1
plant_uptake: 10
plant_growth: 0.05
plant_senescence: 0.1
plant_uptake_saturation: 3
water_infilt_saturation: 5
veg_mass_per_cell: 90
fraction_plant_cells: 0.9<jupyter_text>Setting the initial conditionsFor the starting pattern of vegetation, we can either use a pre-defined pattern (e.g. taken from a satellite image), or define a random pattern, where the fraction of grid cells with non-zero plant biomass is defined in the configuration file (default is 10%).If we use a random pattern though, we first need to set the daily rainfall, in order to get values for the equilibrium quantity of soil water in the bare-soil cells.<jupyter_code>pg.set_rainfall(1.4)
pg.set_random_starting_pattern()<jupyter_output><empty_output><jupyter_text>Now we use this pattern to set the initial values for soil and surface water, and plant biomass for all cells:<jupyter_code>pg.initial_conditions()
pg.plot_image()<jupyter_output><empty_output><jupyter_text>Running the simulationNow we can evolve the simulation by 1000 steps.  In each step, some quantity of surface water will enter the soil, some quantity of soil water will be uptaken by the plants, there will be some diffusion from cell-to-cell, and the plants will gain or lose biomass according to the these values and the parameters in the configuration.<jupyter_code>pg.evolve_pattern(1000)<jupyter_output>Doing 1000 steps with rainfall 1.4mm<jupyter_text>Let's see what the plant biomass looks like now:<jupyter_code>pg.plot_image()<jupyter_output><empty_output><jupyter_text>Quantifying the connectedness of the patterned vegetationIn order to use the subgraph centrality method, we need to convert this into a binary array.  To do this, we define a threshold of (max_val + min_val)/2 and set everything above this to be 0 (i.e. a black pixel, which we define as vegetation), and everything below it as 255 (i.e. a white pixel, in 8-bit color).<jupyter_code>binary_array = pg.make_binary()
binary_array<jupyter_output><empty_output><jupyter_text>The network centrality code is implemented in the "subgraph_centrality.py" module in pyveg.<jupyter_code>from pyveg.src.subgraph_centrality import *
feature_vec, _ = subgraph_centrality(binary_array)<jupyter_output><empty_output><jupyter_text>The "feature_vector" is the Euler Characteristic values for the different quantiles of vegetation-covered pixels, ordered by subgraph centrality.<jupyter_code>plt.plot(feature_vec, "ro")<jupyter_output><empty_output><jupyter_text>The single number that we are using to quantify the connectedness is "offset50", which is essentially the slope of the second half of this feature vector.<jupyter_code>offset50 = feature_vector_metrics(feature_vec)["offset50"]
print(offset50)<jupyter_output>-378.0<jupyter_text>Generate a precipitation time-seriesLet's assume we're in a semi-arid landscape, where we get very little rainfall for most months, then between 1 and 2mm per day in June, July, August<jupyter_code>n_months = 36
months = np.array(list(range(n_months)))
daily_rainfall = np.zeros(n_months)
for i in range(len(daily_rainfall)):
    if i%12 == 5: # May
        daily_rainfall[i] = 1.0
    elif i%12 == 6: # June
        daily_rainfall[i] = 1.3
    elif i%12 == 7: # July
        daily_rainfall[i] = 1.7
    elif i%12 == 8: # August
        daily_rainfall[i] = 1.4
    else:
        daily_rainfall[i] = 0.8

    
plt.plot(months,daily_rainfall)<jupyter_output><empty_output><jupyter_text>Set a starting vegetation patternThis time let's load a starting pattern from a csv file - this is one of the "standard" labyrinth patterns we've been using in testing:<jupyter_code>pg.set_starting_pattern_from_file(os.path.join("..","..","pyveg"
                                               ,"testdata",
                                               "binary_labyrinths_50.csv"))
pg.plot_image()<jupyter_output>/Users/nbarlow/anaconda3/envs/ds4sdenv/lib/python3.7/site-packages/pyveg/src/pattern_generation.py:392: UserWarning: Matplotlib is currently using ps, which is a non-GUI backend, so cannot show the figure.
  plt.show()<jupyter_text>Run the simulationNow we can step through this 3-year period, setting the rainfall according to the time series we generated above, and seeing how the vegetation responds:<jupyter_code>offset50_vals = np.zeros(n_months)
for i in range(n_months):
    # set the rainfall for month i
    pg.set_rainfall(daily_rainfall[i])
    # evolve the pattern for 30 days
    pg.evolve_pattern(30)
    # get a binary array representing the vegetation pattern
    binary_array = pg.make_binary()
    # calculate the subgraph centrality, and set "offset50"
    feature_vec, _ = subgraph_centrality(binary_array)
    offset50 = feature_vector_metrics(feature_vec)["offset50"]
    print("offset50 value is: {}".format(offset50))
    offset50_vals[i] = offset50<jupyter_output>Doing 30 steps with rainfall 0.8mm
offset50 value is: -451.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -427.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -508.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -802.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -1109.0
Doing 30 steps with rainfall 1.0mm
offset50 value is: -1115.0
Doing 30 steps with rainfall 1.3mm
offset50 value is: -1046.0
Doing 30 steps with rainfall 1.7mm
offset50 value is: -995.0
Doing 30 steps with rainfall 1.4mm
offset50 value is: -1015.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -1033.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -1007.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -946.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -807.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -756.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -755.0
Doing 30 steps with rainfall 0.8mm
offset50 value is: -763.0
Doing 30 steps wit[...]<jupyter_text>Results<jupyter_code>plt.plot(offset50_vals)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07929454242466324, 'fraction_numerical': 0.04568809887515623, 'mean_word_length': 5.031825795644891, 'pattern_counts': {'":': 0, '<': 56, '<?xml version=': 0, '>': 56, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '15', 'dir': 'jupyter-structured-clean-dedup', 'id': '33103', 'license': 'MIT', 'n_tokens_mistral': 2541, 'n_tokens_neox': 2226, 'n_words': 792, 'path': 'notebooks/tutorials/tutorial_simulate_patterned_vegetation.ipynb', 'repo_name': 'usethedata/monitoring-ecosystem-resilience', 'size': '7201'}