[R] trying to obtain same nls parameters as in example

[Pedro Mardones]

Dear R-users;

I'm working with a a dataset that was previously used to fit a
nonlinear model of the form:

Y ~ a * (1 + b * log(1 - c * X^d))

The parameters published elsewhere are:

a = 1.758863, b = .217217, c = .99031, and d = .054589

However, there is no way I can replicate this result. I've tried
several options (including SAS) w/o success.

The data is:


X <- c(0.0048,0.0145,0.0338,0.0628,0.1353,0.2077,0.2802,0.3527,0.4251,0.4976,0.57,0.6425,0.715,0.7874,0.8599,0.8841,0.0049,0.0148,0.0345,0.064,0.1379,0.2118,0.2857,0.3596,0.4335,0.5074,0.5813,0.6552,0.7291,0.803,0.8374,0.0052,0.0156,0.0365,0.0677,0.1458,0.224,0.3021,0.3802,0.4583,0.5365,0.6146,0.6927,0.7708,0.8333,0.0049,0.0148,0.0345,0.064,0.1379,0.2118,0.2857,0.3596,0.4335,0.5074,0.5813,0.6552,0.7291,0.803,0.8621,0.0055,0.0166,0.0387,0.0718,0.1547,0.2376,0.3204,0.4033,0.4862,0.5691,0.6519,0.7348,0.8177,0.8564,0.0052,0.0155,0.0363,0.0674,0.1451,0.2228,0.3005,0.3782,0.456,0.5337,0.6114,0.6891,0.7668,0.8446,0.886,0.0057,0.0172,0.0402,0.0747,0.1609,0.2471,0.3333,0.4195,0.5057,0.592,0.6782,0.7644,0.8506,0.0055,0.0164,0.0383,0.071,0.153,0.235,0.3169,0.3989,0.4809,0.5628,0.6448,0.7268,0.8087,0.8525,0.0055,0.0166,0.0387,0.0718,0.1547,0.2376,0.3204,0.4033,0.4862,0.5691,0.6519,0.7348,0.8177,0.8729,0.0058,0.0174,0.0407,0.0756,0.1628,0.25,0.3372,0.4244,0.5116,0.5988,0.686,0.7733,0.8605,0.0056,0.0168,0.0391,0.0726,0.1564,0.2402,0.324,0.4078,0.4916,0.5754,0.6592,0.743,0.8268,0.8659,0.0057,0.017,0.0398,0.0739,0.1591,0.2443,0.3295,0.4148,0.5,0.5852,0.6705,0.7557,0.8409,0.8636)

Y <- c(1.2948,1.2032,1.1355,1.008,0.9044,0.8127,0.7928,0.7171,0.6454,0.5896,0.5299,0.5219,0.4064,0.2988,0.1793,0.1474,1.2963,1.1806,1.088,1,0.8796,0.8241,0.7685,0.6991,0.6204,0.5833,0.5509,0.4306,0.375,0.2778,0.1898,1.1818,1.1483,1.0622,0.9928,0.8995,0.8134,0.7416,0.6746,0.5742,0.6507,0.5455,0.4258,0.3493,0.2536,1.3374,1.1399,1.0782,0.9877,0.8519,0.786,0.7572,0.716,0.6132,0.5844,0.5432,0.3992,0.358,0.2469,0.1811,1.2687,1.1592,1.0697,1,0.9154,0.8458,0.806,0.7562,0.6468,0.5721,0.5572,0.3682,0.2687,0.194,1.2529,1.1518,1.0389,0.9922,0.8911,0.821,0.751,0.7237,0.6109,0.537,0.5175,0.4163,0.3113,0.1946,0.1518,1.1588,1.0773,1.0644,0.9828,0.8455,0.7811,0.7167,0.6953,0.5751,0.515,0.4464,0.2876,0.1931,1.2188,1.1339,1.1205,1.0045,0.8884,0.7946,0.7813,0.7188,0.683,0.6205,0.6161,0.3482,0.2411,0.1741,1.2845,1.1255,1.0418,1.0084,0.8912,0.7741,0.7406,0.6611,0.6025,0.5397,0.4895,0.3473,0.2469,0.1632,1.425,1.155,1.065,0.99,0.88,0.77,0.725,0.635,0.545,0.515,0.375,0.285,0.21,1.2396,1.2031,1.0885,1.0052,0.875,0.7969,0.7448,0.6823,0.6198,0.5781,0.4583,0.3073,0.2396,0.224,1.2473,1.2151,1.0591,1.0108,0.914,0.7903,0.7312,0.6613,0.5914,0.5323,0.4462,0.3226,0.2688,0.2366)

test <- data.frame(X = X, Y = Y)

and I tried

fit0 <- nls(Y ~ a * (1 + b * log(1 - c * X^d)), data = test, start =
list(a = 1.75, b = .22, c = .99, d = .005))

and

library(minpack.lm)
fit2 <- nlsLM(Y ~ a * (1 + b * log(1 - c * X^d)), data = test, start =
list(a = 1.75, b = .22, c = .99, d = .005), control =
nls.lm.control(maxiter = 100))

I'm not an expert in fitting models so I'm wondering if there is
something I'm missing. I would appreciate your comments/hints.

Thanks

PM