%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Matlab code to simulate seed dispersal kernels from GPS data of % foraging mallards. % % Supplement to the 2017 Journal of Ecology article: % 'Seed dispersal distributions resulting from landscape-dependent daily % movement behaviour of a key vector species'. % Authors: Erik Kleyheeg, Jelle Treep, Monique de Jager, Bart A. Nolet, % Merel B. Soons % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Initialization %%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear; clc; close all %% Data files % sites: 1. Oud Alblas, 2. Enterveen, 3. Terra Nova, 4. Juliusput site = 1; % read GPS data from CSV file if site == 1 % Oud Alblas GPS_init = readtable('Mallard_GPS_tracks_9days.csv'); GPS = GPS_init(3270:end,:); % Select rows roost = 27; % Identifier of roosting site field = 'Oud Alblas'; % Name of the site used in figures elseif site == 2 % Enterveen GPS_init = readtable('Mallard_GPS_tracks_9days.csv'); GPS = GPS_init(1:3268,:); % Select rows roost = 11; % Identifier of roosting site field = 'Enterveen'; % Name of the site used in figures elseif site == 3 % Terra Nova GPS_init = readtable('Mallard_GPS_tracks_9days_JP_TN.csv'); GPS = GPS_init(1:4021,:); % Select rows roost = 18; % Identifier of roosting site field = 'Terra Nova'; % Name of the site used in figures elseif site == 4 % Juliusput GPS_init = readtable('Mallard_GPS_tracks_9days_JP_TN.csv'); GPS = GPS_init(4022:end,:); % Select rows roost = 6; % Identifier of roosting site field = 'Juliusput'; % Name of the site used in figures end % specify filenames for figures and output data filenamepdf = strcat('Kernels_',num2str(site),'.pdf'); % PDF file for printing figures filename = strcat('Kernels_',num2str(site),'.mat'); % .mat file for saving output data %% Model parameters Seed_size = 2; % Volume 0.2, 2 or 20 mm3 Theta = 1/(0.56456-0.0428 * log(Seed_size)); % Scale parameter of gamma distributed retention times k = 2.7; % Shape parameter of gamma distributed retention times Max_RT = 48; % Maximum retention time [hours] Survival=-0.039 * log(Seed_size) + 0.2121; % Seed size dependent survival parameter %% Initialize variables % Read variables from GPS data table ID = table2array(GPS(:,1)); % GPS ID date_vect = datetime(table2array(GPS(:,5)),'InputFormat','dd-M-yyyy'); % Create date vector time_vect = datetime(table2array(GPS(:,6)),'InputFormat','HH:mm:ss'); % Create time vector Year = year(date_vect); Month = month(date_vect); Day = day(date_vect); % Date identifiers of GPS fixes Hour = hour(time_vect); Minute = minute(time_vect); Second = zeros(length(time_vect),1); % Time identifiers of GPS fixes Date_Time = datetime(Year,Month,Day,Hour,Minute,Second); % Translate GPS time-fix to date-time format Int_dur = Date_Time(2:length(Date_Time))-Date_Time(1:(length(Date_Time)-1)); % GPS interval [minutes] xy = table2array([GPS(:,7) GPS(:,8)]); % Coordinates (lat lon) core = table2array(GPS(:,20)); % Core areas where the mallard tod = table2array(GPS(:,15)); % Initialize foraging variables NumSeeds = zeros(length(Date_Time),1); % Number of (surviving) seeds inside the mallard each timestep Intake = zeros(length(Date_Time),1); % Seed intake of the mallard each timestep Excretion = zeros(length(Date_Time),1); % Seed excretion of the mallard each timestep Int = minutes(Int_dur); % Convert duration to real numbers [unit minutes] Dist_prob = zeros(100,1); % Initializes dispersal kernel %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Dynamic part %%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ii = 1; % Iterator that resets to zero when simulation of new individual starts tic for i = 1:length(xy(:,1))-1 %% Seed part % Seed intake if ((length(tod{i}) ~= 3) && core(i) ~= roost) % If it is nighttime and mallard is not at its roost Intake(i) = 1000*Int(i)*Survival/60; % Intake rate [# seeds per time step] end % Seed excretion rate if ii > 1 if sum(ID{i} == ID{i+1}) == 3 && (site ~= 3 || i ~= 452) % When next fix is new individual or when site = 3 and i = 452 (time gap in database) % Retention curve jjj = min(ii-1,Max_RT*4+2); % Specifies length of the excretion vector (*4 because most timesteps are 15 minutes; +2 because some are 13-14 minutes) cum_hr = cumsum(Int(i-1:-1:(i-jjj))/60); % Cumulative time from all GPS timesteps within preliminary retention period (may be longer due to long GPS intervals) abc = find(cum_hr < 48); % Identifier to select only GPS timesteps within retention maximum retention time Retention = gampdf(cum_hr(abc),k,Theta); % Makes retention curve % Calculate excretion Excretion(i) = sum(Retention.*Intake(i-abc))*Int(i)/60; % Excretion rate [# seeds per time step] Excretion_vec = Retention.*Intake(i-abc).*Int(i)/60; % Calculates the excretion at this timestep from seed eaten at previous timesteps seperately % Update Number of seeds in mallard NumSeeds(i) = NumSeeds(i-1) + Intake(i) - Excretion(i); % Updates state [# seeds] %% Movement part jjj=min(i,length(Excretion_vec)); % Specifies length of the excretion vector in case the excretion vector has fewer cells than Max_RT*4+2 Distance_vec = dist_vec_ED(xy(i+1,2),xy(i+1,1),... xy(i:-1:(i+1-jjj),1),xy(i:-1:(i+1-jjj),2),jjj); % Calculates distances to previous relevant locations with function dist_vec_ED % Distance probability (Dist_prob as a function of Dist_class Dist_class=linspace(100,10000,100); % Creates distance bins of 100 meter for dispersal kernel for iii=1:jjj % match distance jj=1; while Dist_class(jj) < Distance_vec(iii) && jj~=100 jj=jj+1; end % update Dist_classes Dist_prob(jj) = Dist_prob(jj)+Excretion_vec(iii); % Puts excreted seeds in the distance bins of dispersal kernel end else %if sum(ID{i} == ID{i+1}) ~= 3 ii=0; end else NumSeeds(i) = 0; end ii=ii+1; end toc %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Visualization %%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% figure plot(Dist_class,Dist_prob/sum(Dist_prob),'linewidth',2) axis([0 10000 0 0.3]) xlabel('Distance [m]','fontweight','bold','fontsize',15) ylabel('Probability','fontweight','bold','fontsize',15) title(strcat('\bf',field, {' - Seed Volume = '},num2str(Seed_size), ' [mm^3]')) h=gcf; set(h,'PaperPositionMode','auto'); set(h,'PaperOrientation','landscape'); set(h,'Position',[50 50 500 400]); print(gcf, '-dpdf', filenamepdf) save(filename)