This is the code for the separately fitting.
G=2.529e-4;
R=[0.1:0.1:0.4 0.6:0.2:1.8 2.1 2.4:0.4:3.2 3.8 4.8 6 8.4 12.8 20 36];
pks_locs1=[136 123.5 121 120.5 124 119 118 114 111.5 109.5 108 103 103.5 101.5 98 97 96 90 84 77.5 74.5 56 63.5];
pks_locs2=[903 920 951.5 964.5 944 1006 1026 1050.5 1105.5 1112.5 1122.5 1169 1177 1201 1229 1247 1235.5 1321 1354 1404.5 1444 1499.5 1474.5];
init=[100 1000];
fitted1=lsqcurvefit(@(param,x) delta(param,x,pks_locs1,G),init,R,pks_locs1(2:end));
fitted2=lsqcurvefit(@(param,x) delta(param,x,pks_locs2,G),init,R,pks_locs2(2:end));
figure(1)
subplot(1,2,1),plot([0 R],pks_locs1,'ko',R,delta(fitted1,R,pks_locs1,G),'b-')
str1=['K_a H8 = ',num2str(round(fitted1(2)))];
text(18,(min(pks_locs1)+max(pks_locs1))/2,str1)
legend('Data','Fit','Location','northeast')
title('Fit for H8 chemical shifts')
subplot(1,2,2),plot([0 R],pks_locs2,'ko',R,delta(fitted2,R,pks_locs2,G),'b-')
str2=['K_a H3 = ',num2str(round(fitted2(2)))];
text(18,(min(pks_locs2)+max(pks_locs2))/2,str2)
legend('Data','Fit','Location','northeast')
title('Fit for H5 chemical shifts')