For context, I have calculated 3 functional diversity indices from the potentially non-affected fraction (PNAF) of species at different dissolved oxygen concentrations. In order to fit a regression I am using the nls function with the formula provided by the peper which published the data. When reproducing their results, the fit I have used worked well with the following code;
globalfit <- nls(PNAF ~ 1/(1 + exp(-(log10(LOEC)-a)/b)),
data = global,
start = c(a=0.3 , b=0.2))
I calculated the PNAF of each functional indice using a code similar to the following
#Functional Divergence of fish
FDiv_fish <- data.frame(Divergence= FD_calc_fish$FDiv, LOEC= AllLOEC_fish)
#Calculate the fraction of functional evenness compared to the functional richness at 100% population
FDiv_fish$Divergence_frac[1] <- FDiv_fish$Divergence[1]/max(FDiv_fish$Divergence)
for (i in 2:nrow(FDiv_fish)) {
FDiv_fish$Divergence_frac[i] <- FDiv_fish$Divergence[i]/FDiv_fish$Divergence[nrow(FDiv_fish)]
}
#Remove empty rows
FDiv_fish <- na.omit(FDiv_fish)
# Calculate species richness
SRic_fish <- rowSums(pnaf_matrix_fish)
SRic_fish <- data.frame(LOEC=AllLOEC_fish, SRic = SRic_fish)
SRic_fish$SRic_frac <- SRic_fish$SRic/max(SRic_fish$SRic)
The resulting dataframe for the Functional Divergence indice is the following
FDiv_fish = {Divergence = c(0.889950737228376, 0.955047157345304,
0.952619624560615, 0.935368469313814, 0.940003820276155, 0.942700206928084,
0.942650807180955, 0.940289064256148, 0.943752991360123, 0.941559360809695,
0.924050406506512, 0.926648389394962, 0.876447636275341, 0.881548684676512,
0.886245224649011, 0.88218834971409, 0.885314645756991, 0.881492706719018,
0.876662397962403, 0.886349294318137, 0.883035022314938, 0.8848327450232,
0.881672831746863, 0.884156297551429, 0.881140906361141, 0.883534593605944,
0.885833771865729, 0.88804392997383, 0.891027945817535, 0.888186735006869,
0.893358872293167, 0.886332319236113, 0.888904363556696, 0.887038353743163,
0.884779304420887, 0.897837333190895, 0.894509967650634, 0.891991192039641,
0.899398165601533, 0.89549768133025, 0.885348372920148, 0.883072158209693,
0.880708851749932, 0.870967680069655, 0.872686632273436, 0.870231057490782,
0.871863360688521, 0.873455110028267, 0.888942354125257, 0.882187190419575,
0.88359490643209, 0.832093698279846, 0.833758685259681, 0.837055666849556,
0.838621975644644, 0.84572076038212, 0.853563980765286, 0.852846371924774,
0.848690020275423, 0.848019579107283, 0.897076657728166, 0.886592898915272,
0.875127667801321, 0.872745593426437, 0.870399511052835, 0.873366955360124,
0.864058267778409, 0.852842506107893, 0.854601294578928, 0.852417777332101,
0.854148580874727, 0.852008859522321, 0.854982653095199, 0.853592877986531,
0.852255836650848, 0.844261486012797, 0.842393713964801, 0.851535546465613,
0.845919707924595, 0.83799275396368, 0.839680366487134, 0.848107377439518,
0.850662066070129, 0.843754962494148, 0.840898540618949, 0.839502768752552,
0.841098849929625, 0.836293864179716, 0.831588394725384, 0.828942721595275,
0.831778862465456, 0.827309791131041, 0.822929138401967, 0.82559644603514,
0.823170938369431, 0.820778350271261, 0.818418016069642, 0.816089287920916,
0.818746480527328, 0.81984077135057, 0.817565361142596, 0.819146561499075,
0.820075084732484, 0.820993582920095, 0.820910837777695, 0.820524854216395,
0.820811742816053, 0.822601708617993, 0.822873768873091, 0.823752214457896,
0.82401613553359, 0.824880160298616, 0.825735443483524, 0.825986203258117,
0.826234629738616, 0.826480755340803, 0.827312025301094, 0.827550864133838,
0.832160021463369, 0.822561048178852, 0.81328819769379), LOEC = c(1.67053607,
1.75484884, 1.80915198, 1.83058743, 1.86345512, 1.86774221, 1.8720293,
1.89060669, 1.89489378, 1.90775505, 1.9149002, 1.95634207, 1.96920334,
2.000642, 2.00778715, 2.02636454, 2.07780962, 2.0863838, 2.09495798,
2.09924507, 2.11210634, 2.13068373, 2.16640948, 2.19070299, 2.21070941,
2.21213844, 2.22785777, 2.23357389, 2.26929964, 2.27501576, 2.27644479,
2.3578995, 2.36075756, 2.36790271, 2.39362525, 2.51366377, 2.51937989,
2.69086349, 2.75945693, 2.83233746, 2.85806, 2.8723503, 2.91236314,
2.93379859, 2.93951471, 3.08384674, 3.13957891, 3.35107535, 3.40394946,
3.43967521, 3.44539133, 3.50826865, 3.59829754, 3.63973941, 3.65117165,
3.74834569, 3.95698407, 4.00414206, 4.03558072, 4.03986781, 4.1298967,
4.21992559, 4.24136104, 4.2727997, 4.39426725, 4.40712852, 4.50716062,
4.52287995, 4.58004115, 4.58861533, 4.71008288, 4.71436997, 4.75724087,
4.77010214, 4.8015408, 4.81011498, 4.83440849, 4.84155364, 4.85155685,
4.95158895, 5.10449516, 5.14022091, 5.14736606, 5.16451442, 5.25311428,
5.28598197, 5.29312712, 5.38315601, 5.67610716, 5.71612, 5.72326515,
5.77471023, 5.78900053, 6.00764212, 6.02621951, 6.08623877, 6.13196773,
6.14768706, 6.16197736, 6.18341281, 6.20913535, 6.41062858, 6.58068315,
6.59354442, 6.68071525, 6.70786682, 6.72215712, 6.85076982, 7.02511148,
7.04368887, 7.45810757, 7.5166978, 7.55099452, 7.94254874, 8.00113897,
8.02114539, 8.07116144, 8.1026001, 8.56274776, 9.52162689, 10.24042898
), Divergence_frac = c(1.09426245179996, 1.17430347575865, 1.17131863865961,
1.1501070247499, 1.15580654304549, 1.15912195652324, 1.15906121575844,
1.15615727232055, 1.16041643544845, 1.15771919902396, 1.13619060147043,
1.13938501999983, 1.0776593571143, 1.08393148600494, 1.0897062408653,
1.08471800305928, 1.0885620229919, 1.08386265682772, 1.07792342302313,
1.0898342024777, 1.08575905173458, 1.08796948920725, 1.0840841343167,
1.08713774533873, 1.08343009139904, 1.08637331281992, 1.08920032822025,
1.09191788654012, 1.09558696209313, 1.09209347624307, 1.09845301435141,
1.08981333031692, 1.09297585539459, 1.0906814537067, 1.08790377990216,
1.10395962432119, 1.0998683740735, 1.09677134694568, 1.10587878706702,
1.10108284353514, 1.08860349311683, 1.08580471315554, 1.08289884723192,
1.07092133211748, 1.07303491523433, 1.07001559835549, 1.07202263989669,
1.07397981736989, 1.09302256770232, 1.08471657761807, 1.08644747204947,
1.02312280030545, 1.02517002905481, 1.02922391991321, 1.0311498162923,
1.03987831469865, 1.0495221536298, 1.04863979871238, 1.04352924668281,
1.04270488802368, 1.10302431570011, 1.09013373294897, 1.07603635498817,
1.07310741247844, 1.07022272488522, 1.07387142446761, 1.06242568160781,
1.04863504539506, 1.05079761024725, 1.04811280890251, 1.05024096414629,
1.04761001320114, 1.05126651969085, 1.04955768497198, 1.04791369045752,
1.03808402532686, 1.035787456837, 1.04702803862183, 1.04012293590801,
1.03037614014312, 1.03245118872767, 1.04281284278373, 1.04595402771406,
1.03746121594627, 1.03394902692976, 1.03223281873891, 1.0341953225372,
1.02828722530483, 1.02250149096407, 1.01924843363752, 1.02273568560825,
1.01724062082422, 1.01185427347343, 1.01513393207506, 1.01215158501459,
1.00920971507851, 1.00630750377344, 1.00344415452612, 1.00671137592924,
1.00805689013484, 1.00525909937084, 1.00720330606284, 1.00834499634685,
1.00947436007083, 1.0093726185939, 1.00889802230393, 1.00925077376457,
1.01145167352805, 1.0117861924057, 1.01286630839323, 1.01319081952772,
1.01425320401513, 1.01530484006165, 1.01561316837049, 1.0159186277159,
1.01622125795557, 1.01724336790706, 1.01753703850676, 1.02320434972878,
1.01140167841038, 1), category = c("LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51", "LOEC<3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51",
"LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51", "LOEC>=3.51")))}
When trying to apply the formula used for the SSD curves with the results from the functional diversity indices, I have my third indice which shows an increase in divergence (when dissolved oxygen contrations diminish) which in this case transaltes to an decreasing trend. I cannot seem to find a way to make the nls modelling work for this one.
Here is the code I have tried so far
# PNAF of Functional richness
nls_FRic_fish <- nls(Richness_frac ~ 1/(1 + exp(-(log10(LOEC)-a)/b)),
data = FRic_fish,
start = c(a=0.3 , b=0.1))
# PNAF of Functional Evenness
nls_FEve_fish <- nls(Evenness_frac ~ 1/(1 + exp(-(log10(LOEC)-a)/b)),
data = filter(FEve_fish, LOEC>=3),
start = c(a=0.3 , b=0.1),
trace=TRUE)
# PNAF of Functional Divergence
nls_FDiv_fish <- nls(Divergence_frac ~ log(1/(1 + exp(-(log10(LOEC)-a)/b))),
data = FDiv_fish,
start = c(a=0.3 , b=0.1),
trace = TRUE)
# Create a new variable to categorize LOEC values (identified 2 main behaviours)
FEve_fish$category <- ifelse(FEve_fish$LOEC < 3, "LOEC<3", "LOEC>=3")
FDiv_fish$category <- ifelse(FDiv_fish$LOEC < 3.51, "LOEC<3.51", "LOEC>=3.51")
# Plot
tiff('scatter.fish.tiff', width = 12, height = 4, units = 'in', res = 300)
plot1 <- ggplot(FRic_fish, aes(x = LOEC, y = Richness_frac)) +
geom_point(size = 0.8) +
geom_function(fun= function(x) 1/(1 + exp(-(log10(x)-coef(nls_FRic_fish)['a'])/coef(nls_FRic_fish)['b'])),
color="blue",
size=0.4) +
labs(x = expression("Dissolved oxygen concentration " (~ mgO[2] ~ L^{-1})), y = "PNAF of functional richness") +
xlim(0, 10.5) +
theme_bw() +
theme(panel.border = element_blank())
plot2 <- ggplot(FEve_fish, aes(x = LOEC, y = Evenness_frac)) +
geom_point(size = 0.8) +
geom_function(fun= function(x) 1/(1 + exp(-(log10(x)-coef(nls_FEve_fish)['a'])/coef(nls_FEve_fish)['b'])),
color="blue",
size=0.4) +
labs(x = expression("Dissolved oxygen concentration " (~ mgO[2] ~ L^{-1})), y = "PNAF of functional evenness") +
xlim(0,10.5) +
theme_bw() +
theme(panel.border = element_blank())
plot3 <- ggplot(FDiv_fish, aes(x = LOEC, y = Divergence_frac)) +
geom_point(size = 0.8) +
geom_function(fun= function(x) 1/(1 + exp(-(log10(x)-coef(nls_FDiv_fish)['a'])/coef(nls_FDiv_fish)['b'])),
color="blue",
size=0.4) +
labs(x = expression("Dissolved oxygen concentration " (~ mgO[2] ~ L^{-1})), y = "PNAF of functional divergence") +
xlim(0,10.5) +
theme_bw() +
theme(panel.border = element_blank())
grid.arrange(plot1, plot2, plot3, ncol = 3)
dev.off()
When attempting to produce the nls_SDiv_fish model , the usual error pop up such as "step factor 0.000488281 reduced below 'minFactor' of 0.000976562", "singular gradient", "singular gradient matrix at initial parameter estimates", etc.
I've tried using different starting parameters including negative values, I've tried using different algorithms and different versions of nls (nls2 and nlsLM) with no solution. Does anyone know how to troubleshoot this? Tell me if you need to see more of my code.