clear load diseaszdat.mat %2 removes current mammal load grzeronobatsnoland.mat % no size threshold, remove bats, no land mask % calculate precipitatoin and seasonality dataset %https://neo.sci.gsfc.nasa.gov/view.php?datasetId=TRMM_3B43M&year=2014 for i = 1:360 for j=1:180 pricip2z(i,j) = nanmean(precipz(i,j,:)); end end precip1z=pricip2z'; precip3z=[precip1z(:,181:end) precip1z(:,1:180)]; %calculate megafauna variables %change in home range and fecal diff HRchav = ((HR+HRmega)./(N+Nummega))./(HR./N);%past/current home range Dismovchav = ((Dismov+Dismovmega)./(N+Nummega))./(Dismov./N);%past/current fecal diffusivity %past dispersal HRchme = ((HR+HRmega)./(N+Nummega));%past home range Dismovchme = ((Dismov+Dismovmega)./(N+Nummega));%past fecal diffusivity %current dispersal HRch = (HR./N);%current home range Dismovch = (Dismov./N);%current fecal diffusivity Mweights = (Msummega+Msum); %biomass weighted species richness (N*weight) %choose model - 1=HR, 0=FD mdty = 0; for i = 1:180 for j = 1:360 %This produces a EID map with no megafauna extinctions dismap2me(i,j) = -0.11597006 + 0.12548595*log10(popden(i,j)) + 0.00044051 *precip3z(i,j) ; if mdty==1 %model HR from table 2 dismap2(i,j) = -0.11184392 + 0.12645528*log10(popden(i,j)) + 0.00174199*HRchav(i,j) + 0.00042810 *precip3z(i,j) ; %high and low value of HR value jmap(i,j) = ( 0.00174199- 0.00054924 )*HRchav(i,j); jmap2(i,j) = ( 0.00174199+ 0.00054924 )*HRchav(i,j); else %model FD from table 2 dismap2(i,j) = -0.11597006 + 0.12548595*log10(popden(i,j)) + 0.00258366 *Dismovchav(i,j) + 0.00044051 *precip3z(i,j) ; %high and low value of FD value jmap(i,j) = (0.00258366-0.00069726) *Dismovchav(i,j); jmap2(i,j) = (0.00258366+0.00069726) *Dismovchav(i,j); end end end dismap22=dismap2(:); %without JID, the prediction goes to less than zero - we force it back to %zero below md = min(dismap22); dismap2me=dismap2me-md; dismap2=dismap2-md; mapdiff = (dismap2-dismap2me); nanmean(mapdiff(:)); jmapz = landz.*jmap; jmap2z = landz.*jmap2; jmapz(jmapz==0)=NaN; jmap2z(jmap2z==0)=NaN; hivl = (nanmean(jmapz(:)))/nanmean(dismap2(:)); lovl = (nanmean(jmap2z(:)))/nanmean(dismap2(:)); hivlc = hivl/5; lovlc = lovl/5; hivlcd=hivlc*145; lovlcd=lovlc*145; zzz = [hivl lovl; hivlc lovlc; hivlcd lovlcd] % This produces fig 4b impact of extinctions - where we subtract best EID map from map with % megafauna x2=[mapdiff(:,181:end) mapdiff(:,1:180)]; %this produces the best EID map % x2=[dismap2(:,181:end) dismap2(:,1:180)]; geoidR = makerefmat('RasterSize', size(x2), 'Latlim', [-90 90], 'Lonlim', [0 360]); h1 = worldmap('world'); % manually delete long lines land = shaperead('landareas.shp', 'UseGeoCoords', true); %geoshow('landareas.shp', 'FaceColor', [0.5 1.0 0.5]); axis off geoshow(flipud(x2), geoidR, 'DisplayType', 'texturemap'); geoshow([land.Lat],[land.Lon], 'Color', 'red') colorbar('location', 'southoutside'); %title('Bankrupt Countries ', 'FontSize', 10); hcb = colorbar('horiz'); set(get(hcb,'Xlabel'),'String','mol m-2 yr-1') colormapeditor