From: [Permission pending] Newsgroups: sci.math Subject: Help: Curve fitting to data Date: 30 Jan 1995 11:26:15 GMT Hi I could do with some advice on fitting curves to experimental data. I've got data that looks like this x y1 y2 y3 2.05 8.51 1.98 2.36 2.23 8.51 1.96 2.34 2.53 8.54 1.96 2.31 2.83 8.64 1.96 2.28 3.23 9.21 2.03 2.24 3.94 9.21 2.03 2.24 4.91 9.54 2.06 5.91 9.63 2.06 2.12 and I want to fit the curve y = (a + b*10^(x - c)) __________________ (1 + 10^(x - c)) to it. I can fit the curve to each set of data (y1, y2, y3) indidually no problem. However, each set of data should have the same value of c, so what I want is to be able to fit all the data to obtain a single value of c. I've normalised and fitted the data, but you get bad scatter, which is a result of the normalisation (the normalisation emphaises experimental scatter, and makes it look awful). What I'dlike to do would be to fit all three sets of data simultaneously, to get a value of (a,b) for each of y1, y2 and y3, and a value of c obtained for the three sets as a whole. Is this possible, and if it is, how do I do it ? (is there a public domain program that can do this ?). Thanks a lot for any advice, [sig deleted -- djr] ============================================================================== Date: Sat, 4 Feb 95 04:22:21 CST From: rusin (Dave Rusin) To: [Permission pending] Subject: Re: Help: Curve fitting to data Newsgroups: sci.math How large is your data set and how definite do you think the 'right' value of c is? If the answers are 'small' and 'very', respectively, this approach might work: As I understand it, you've got a number of 4-tuples (X=10^x, y1, y2, y3) and you seek 7 constants C = 10^(-c), A1 = a1-b1, A2, A3, b1, b2, b3 such that when we look at the functions f1(X) = b1 + A1/(1+C X), f2, and f3, we have y1(i) roughly equal to f1(X(i)) for all i (that is, |y1(i) - f1(X(i))| is small), and likewise for f2 and f3. You have to decide what your measure of aggregate error is; most people go for "small" in the least-squares sense, which would mean we want to choose C, A1, A2, A3, b1, b2, b3 to minimize Sum ( y1(i) - f1(X(i)) )^2 + Sum ( y2-f2)^2 + Sum (y3-f3)^2 (pardon my abbreviation).This hits a minimum when all 7 partials are zero. Taking the partial derivatives with respect to all variables except C gives 3 pairs of linear equations in Ai and bi, which are easily solved (in terms of C). Using these solutions, we can express the last condition (partial deriv. of aggregate error with respect to C) as a function F(C) of C only; we seek C which makes F(C) = 0. Unfortunately this function is terrible -- it involves a lot of terms 1/(1+C X(i) ) added together, and so on; it's a rational function of C with huge degrees in numerators and denominator. However, it's relatively straightforward to express the derivative of F as a similar function of C. Therefore, _if_ the dataset is not too large, we can calculate F(C) and F'(C) for any specific value of C. This lets us solve F(C)=0 by newton's method. Of course, you'd like to know where to start with Newton's method and how rapidly it will converge; this is why I hoped you already had a reasonable sense of the value of C. If it looks like this approach would work for you, let me know and I can help with the details. dave ============================================================================== Date: Mon, 6 Feb 95 17:04:25 CST From: rusin (Dave Rusin) To: [Permission pending] Subject: Re: Help: Curve fitting to data Newsgroups: sci.math In article you write: >On Sat, 4 Feb 1995, Dave Rusin wrote: > >> How large is your data set and how definite do you think the 'right' >> value of c is? If the answers are 'small' and 'very', respectively, this >> approach might work: > >OK - got 3 - 4 curves in a data set, each containing 8 points. I know the >value of c to within 10%, probably within 5%. Is this close enough ? 8 points is nice and small; I'm hoping the value of c is close enough. To recap: transform your datasets replacing x with X=10^x. Then you have datasets like (X(1), y1(1), y2(1), y3(1) ) ... (X(8), y1(8), y2(8), y3(8) ) You wanted y1 to be roughly f1(X) = (a + bCX)/(1+CX) where I have written C for 10^(-c); this is equal to b + (a-b)/(1+CX), which I'll write f1(X) = B + A/(1+CX). For any given C, what A and B are best? Well, you try to minimize Sum( (y1(i) - f1(X(i)))^2, i=1 to 8). The minimum occurs when the partial derivatives are equal to zero. This gives two linear equations in A and B whose solution is to take B and A to be the entries in the vector (1/D) M V, where D=8Sum(w(i)^2) - [Sum(w(i))]^2, and M and V are this matrix and vector: ( Sum(w(i)^2) -Sum(w(i)) )( Sum(y1(i)) ) ( )( ) ( -Sum(w(i)) 8 )( Sum(y1(i)w(i) ) The sums are for i=1 to 8 and I've written w(i) for 1/[1+CX(i)]. It is important to realize that A and B thus chosen depend on the choice of C. Indeed, differentiating wrt C the equations which define A and B gives equations defining A'=dA/dC and B'=dB/dC. The result is that these derivatives are the entries in (1/D) M W where D and M are as above and W is the vector ( Sum( X(i)w(i) ) A ) ( ) ( -Sum( y(i)X(i)w(i)^2) + Sum(X(i) w(i)^2) B + 2 Sum(X(i)w(i)^3) A ) We will need below to take the derivative of f1(X) with respect to C . Recall again that f1(X) = B + A/(1+CX), where B and A themselves depend on C. Thus df1/dC =B' + A'/(1+CX) -AX/(1+CX)^2 . Keeping notation as above, this is B' + A'w(i) -AX(i)w(i)^2 when X=X(i). Now the question is which C is best. I'll assume 'best' means you want to minimize the sum of the squared errors as before, now adding together the (squared) errors for y1, y2, and y3. Again you accomplish this by setting the derivative wrt C equal to zero. Thus C is determined by solving one awful equation: 0=Sum (y1(i) - f1(X(i))).(df1(X(i))/dC) + (same for y2) + (same for y3). Each of the three terms has this form: B'Sum(y(i)) - BB' - AB'Sum(w(i)) + +A'Sum(y(i)w(i)) - A'B Sum(w(i)) - AA' Sum(w(i)^2) -A Sum(X(i)y(i)w(i)^2) + AB Sum (X(i)w(i)^2) + A^2 Sum(w(i)^3X(i)) I will write the equation deterimining C as 0=F(C). My suggestion is to use Newton's method to approximate C to any desired degree of accuracy. I can't _guarantee_ convergence but if indeed you have c 'almost right' you should be lucky enough. (Newton's method roughly doubles the number of decimal places of accuracy with each iteration, once you are close enough to the right answer.) Begin with an estimate C(0) of C. Then keep making more approximations C(n+1) = C(n) - F(C(n))/F'(C(n)). This is no big deal; computationally the big problem is in calculating F(C) and F'(C) for each of these approximations C(n) for C. F'(C) will, like F(C), be a sum of three terms. You need 'only' differentiate the three-line formula above with respect to C. It's not really all that bad. The y(i) and X(i) are constants, and the derivative of w(i) is just -X(i) w(i)^2. Everything else you can get from the product rule except you'll need B'' and A''. Once again these come from differentiating the equations defining B' and A', and once again the solution has the form (B'', A'') = (1/D) M Z with the same M and with a certain vector Z (Ready? Z= 2 N1 (B', A') - N2 (B, A) + N3, where A, B, A', B' are as above, N1 = matrix[0, Sum(X(i)w(i)^2), Sum(X(i)w(i)^2), 2Sum(X(i)w(i)^3)] N2 = matrix[0, 2Sum(X(i)^2w(i)^3), 2Sum(X(i)^2w(i)^3), 6Sum(X(i)^2w(i)^4)] N3 = vector[0, 2Sum(y(i)X(i)^2w(i)^3) ]). So here's the deal. I) Compute the transformed dataset (X(i), y1(i), y2(i), y3(i)) II) Compute your starting guess C = 10^(-c) III) Repeat the following steps as often as you like (perhaps stopping if C has not changed by more than .01%) 1) Compute an additional dataset column w(i)=1/(1+CX(i)) 2) Compute and store additional constants by summing over the dataset. (You'll need this in steps (3)-(8) ) N=Sum(1) (=8) Sum(y(i)) (actually 3 sums: one each for y=y1, y=y2, y=y3) Sum(y(i)w(i)) (ditto) Sum(w(i)) Sum(w(i)^2) Sum(X(i)w(i)^2) Sum(X(i)w(i)^3) Sum(X(i)y(i)w(i)^2) (ditto) Sum(X(i)^2w(i)^3) Sum(X(i)^2w(i)^4) Sum(X(i)^2y(i)w(i)^3) (ditto) [The next 4 steps have to be done once for each 'curve', that is, once with 'y' meaning y1, then with y=y2, etc.] 3) Compute the matrix M and the denominator D as above 4) Compute A and B as above 5) Compute A' and B' 6) Compute A'' and B'' 7) Compute F(C) (the three-line formula). Don't forget to add the resulting values for when y=y1, y=y2, etc. 8) Compute F'(C) (the formula requiring the product rule). Yes, there are many terms. I hope this project merits a long program! 9) Compute a new value of C as C - F(C)/F'(C) 10) loop back (if you decide you need to continue improving C) Again, good luck, and contact me if this looks helpful but insufficient. dave PS1: what's the project? PS2: you shouldn't post private email to a newsgroup. I didn't post in the first place because (1) I don't feel expert enough to know I'm not making a total fool of myself in this particular area and (2) it's not necessary to load up the world's computers with data which is only of interest to you and me. ============================================================================== Date: Tue, 7 Feb 1995 08:33:01 +0000 (GMT) From: [Permission pending] To: Dave Rusin Subject: Re: Help: Curve fitting to data Hi On Mon, 6 Feb 1995, Dave Rusin wrote: > > > >OK - got 3 - 4 curves in a data set, each containing 8 points. I know the > >value of c to within 10%, probably within 5%. Is this close enough ? > > 8 points is nice and small; I'm hoping the value of c is close enough. > > > Again, good luck, and contact me if this looks helpful but insufficient. > This looks excellent - I'll give it a try, see how it goes. I also got sent a mathematica protocol as well, so I'll try both, see how they compare. > dave > > PS1: what's the project? I'm studying proteins in solution with NMR. I did a pH titration to get the pKa's of the carboxylate groups, and each proton acts as a probe for the titration of that residue. I've analysed the data by averaging and normalising the curves, but a simultaneous fitting looked like by far the best way of getting the pKa values. I've got saome more complex models which are for coupled titrations. > PS2: you shouldn't post private email to a newsgroup. I didn't post in > the first place because (1) I don't feel expert enough to know I'm not > making a total fool of myself in this particular area and (2) it's > not necessary to load up the world's computers with data which is only > of interest to you and me. We've just had our version of pine changed; the default is to reply to the newsgroup, and unless you're on your toes (which I'm not always at 8 am) it can escape to the group. Which is very irritating, both for me & for everyone else who has to read it. Thanks again for your help; I'm constantly amazed and impressed by the power of the net to help people across subject boundaries. [sig deletd] ============================================================================== Date: Tue, 7 Feb 95 03:21:44 CST From: rusin (Dave Rusin) To: [Permission pending] Subject: corrected (and shortened) analysis I realized, after writing you this afternoon, that the situation is much simpler than I had imagined it. If you want, you can trash the previous letter since I'll start from scratch here. (Well, I'll recycle parts of that letter). As a bonus: I got it so explicit that I could even write a program to do this (in BASIC, no less!) The good news is that when I try it with data sets which are indeed of the form y=B + A/(1+CX), then it quickly finds the right C if you start reasonably near the right value. The bad news is that when I tried one of your posted data sets -- even after tweaking it a bit -- I couldn't seem to get your data to lie on such a curve. Anyway, I'll send the program separately; here's the derivation, so you can double-check if something seems amiss. .................................... To recap: transform your datasets replacing x with X=10^x. Then you have datasets like (X(1), y1(1), y2(1), y3(1) ) ... (X(8), y1(8), y2(8), y3(8) ) You wanted y1 to be roughly f1(X) = (a + bCX)/(1+CX) where I have written C for 10^(-c); this is equal to b + (a-b)/(1+CX), which I'll write f1(X) = B + A/(1+CX). When convenient, I'll write w(i) for 1/(1+CX(i)). Likewise you wanted y2 to be close to some f2(X), and y3 ~ f3(X). When I need to distinguish these I'll write, say f2(X) = B2 + A2/(1+CX); the C is assumed the same for all three functions. Now, how do you decide when you've got the best C, A1, B1, etc? Usually we try to minimize the 'total error' E = Sum ( (y1(i) - f1(X(i)))^2 ) + similar for y2, y3,...; the sum ranges from 1 to 8 in your case. The nice thing about this E is that to minimize it all you do is set its partial derivatives equal to zero. These partial derivatives are not too bad. If z is any of the variables, we compute the derivative (d/dz) (y1(i) - f1(X(i)))^2 = 2(y1(i)-f(X(i)))( -(d/dz) f1(X(i)) ). If z is any variable except C, A1 or B1, then the last term is zero. If z=B1, the last term is just -1. If z=A1, the last term is -w(i). If z=C, the last term is (-A)[-1/(1+CX(i))^2].X(i) = AX(i)w(i)^2. Writing in f1(X(i)) = B1 + A1 w(i) we see the whole derivative dE/dz is ...for z=B1: 2 Sum ( -y1(i) + B1 + A1 w(i) ) ...for z=A1: 2 Sum ( -y1(i) + B1 + A1 w(i) ) . w(i) (and likewise for z=A2, B2, etc.) Finally, dE/dC is the sum of three parts 2 Sum ( -y1(i) + B1 + A1 w(i) ) .A1 X(i) w(i)^2 and two other similar parts involving A2, B2, y2 and A3, B3, y3. OK, so the conditions that the partials wrt A1 and B1 both vanish is enough to determine A1 and B1 (in terms of C): we must have B1 Sum( 1 ) + A1 Sum( w(i) ) = Sum( y1(i) ) and B1 Sum( w(i) ) + A1 Sum( w(i)^2 ) = Sum( y1(i)w(i) ) , which we can just write in matrix form as M Q1 = V1, where M is a certain 2x2 matrix, Q1=(B1, A1) is the vector we seek, and V1 is the vector of the right-hand sides. This equation is easy to solve. Likewise, we can write and solve two equations determining B2 and A2 (in terms of C) by setting dE/dA2=dE/dB2=0, and so on. There remains one equation to satisfy, namely dE/dC = 0. We have above applied d/dC to something. Extending this we see the derivative becomes the sum of three terms (one each for y1, y2, y3) each of the form -A1 Sum (y1(i) X(i) w(i)^2) + A1 B1 Sum(X(i) w(i)^2) + A1^2 Sum(X(i) w(i)^3) We have already indicated that w(i) and A1 and B1 depend on C, so this equation gives (implicitly, I guess) a condition to be met to make dE/dC=0 What the heck, I'll write it out. We will need to make F(C)=0, where F(C) is the function -A1 Sum (y1(i) X(i) w(i)^2) -A2 Sum (y2(i) X(i) w(i)^2) -A3 Sum (y3(i) X(i) w(i)^2) + (A1 B1 + A2 B2 + A3 B3) Sum(X(i) w(i)^2) + (A1^2 + A2^2 + A3^2) Sum(X(i) w(i)^3) Now, how should we solve the equation F(C)=0? My suggestion is to use Newton's method to approximate C to any desired degree of accuracy. I can't _guarantee_ convergence but if indeed you have c 'almost right' you should be lucky enough. (Newton's method roughly doubles the number of decimal places of accuracy with each iteration, once you are close enough to the right answer.) Begin with an estimate C(0) of C. Then keep making more approximations C(n+1) = C(n) - F(C(n))/F'(C(n)). This is no big deal; computationally the big problem is in calculating F(C) and F'(C) for each of these approximations C(n) for C. F'(C) will, like F(C), be a sum of several terms. We'll have to differentiate the four-line formula above with respect to C. It's not really all that bad. The y(i) and X(i) are constants, and the derivative of w(i) is just -X(i) w(i)^2. Everything else you can get from the product rule except you'll need B' and A' of course. Here's what I get: F'(C) may be written 2 A1 Sum (y1(i) X(i)^2 w(i)^3 ) +2 A2 Sum (y2(i) X(i)^2 w(i)^3) + + 2 A3 Sum (y3(i) X(i)^2 w(i)^3) -A1' Sum (y1(i) X(i) w(i)^2) -A2' Sum (y2(i) X(i) w(i)^2) -A3' Sum (y3(i) X(i) w(i)^2) - 2 (A1 B1 + A2 B2 + A3 B3) Sum(X(i)^2 w(i)^3) + (A1 B1' + A1' B1 + A2 B2' + A2' B2 + A3 B3' + A3' B3) Sum(X(i)w(i)^2) - 3 (A1^2 + A2^2 + A3^2) Sum(X(i)^2 w(i)^4) + 2 (A1 A1' + A2 A2' + A3 A3') Sum(X(i) w(i)^3) Yeah, right, you say. The dependence of A1, etc. on C is hard to see; how will we compute A1' and so on? Actually that's easy. Indeed, differentiating wrt C the equation M Q1 = V1 (the equation which defined A1 and B1 ) gives M Q1' = - M' Q1 + V1' where M is as above and so -M' is the matrix ( 0 Sum( X(i) w(i)^2 ) ) ( Sum( X(i) w(i)^2 ) 2 Sum( X(i) w(i)^3 ) ) and V1' is the column vector (0, -Sum( y(i) X(i) w(i)^2 ) ). ======================================== So here's the deal. I) Compute the transformed dataset (X(i), y1(i), y2(i), y3(i)) II) Compute your starting guess C = 10^(-c) III) Repeat the following steps as often as you like (perhaps stopping if C has not changed by more than 1% or whatever) 1) Compute an additional dataset column w(i)=1/(1+CX(i)) 2) Compute and store additional constants by summing over the dataset. N1=Sum(1) (=8) N2=Sum(w(i)) N3=Sum(w(i)^2) N4,5,6=Sum(y(i)) (actually 3 sums: one each for y=y1, y=y2, y=y3) N7,8,9=Sum(y(i)w(i)) (ditto) N10=Sum(X(i)w(i)^2) N11=Sum(X(i)w(i)^3) N12,13,14=Sum(X(i)y(i)w(i)^2) (ditto) N15=Sum(X(i)^2w(i)^3) N16=Sum(X(i)^2w(i)^4) N17,18,19=Sum(X(i)^2y(i)w(i)^3) (ditto) 3) Compute the matrix M = matrix[ [N1, N2], [N2, N3] ] 4) Solve the vector equations M. [B1, A1] = [N4, N7], M. [B2, A2] = [N5, N8], and M. [B3, A3] = [N6, N9]. 5) Compute the matrix -M' = matrix[ [0, N10], [N10, 2 N11] ] 6) Solve the vector equations M. [B1', A1'] = [ 0, -N12] + (-M')[B1, A1], M. [B2', A2'] = [ 0, -N13] + (-M')[B2, A2], M. [B3', A3'] = [ 0, -N14] + (-M')[B3, A3], 7) Compute F(C) = -A1 N12 - A2 N13 - A3 N14 + + (A1 B1 + A2 B2 + A3 B3) N10 + (A1^2 + A2^2 + A3^2) N11 8) Compute F'(C) = 2 (A1 N17 + A2 N18 + A3 N19) - (A1' N12+ A2' N13+ A3' N14) - 2 (A1 B1 + A2 B2 + A3 B3) N15 + (A1 B1' + A1' B1 + A2 B2' + A2' B2 + A3 B3' + A3' B3) N10 - 3 (A1^2 + A2^2 + A3^2) N16 + 2 (A1 A1' + A2 A2' + A3 A3') N11 9) Compute a new value of C as C(n+1) = C(n) - F(C(n))/F'(C(n)). 10) loop back (if you decide you need to continue improving C) dave ============================================================================== Date: Tue, 7 Feb 95 03:28:28 CST From: rusin (Dave Rusin) To: [Permission pending] Subject: program to find your C Here is a program which seeks to solve the curve-fitting problem you proposed. It runs fine in UBASIC (a BASIC variant for PCs which allows arbitrary precision calculations). Speed of calculation for each iteration (supposedly approximating C) is not a problem, but it's not clear how far from the right C you can start and expect convergence; nor is it clear to what extent you can fit this kind of curve through data which really isn't modelled that way. (Try fiddling with the cooked data to give it data successively further from a true y=B+A/(1+CX). ) The program is pretty lame, but at least it shows I've got all the right signs and so on. You'll surely want to recode it into something prettier (and in a more socially acceptable language) 10 dim Littlex(8),Y1(8),Y2(8),Y3(8),X(8),W(8) 20 'goto 18:' to use actual data 30 'to use cooked data go here: 40 C=1 50 for I=1 to 8 60 X(I)=I:Y1(I)=2+3/(1+C*I):Y2(I)=4+5/(1+C*I):Y3(I)=6+7/(1+C*I) 70 next I 80 goto 140 90 'go here to use actual data 100 for I=1 to 8 110 read Littlex(I),Y1(I),Y2(I),Y3(I) 120 X(I)=10^Littlex(I) 130 next I 140 'OK, one way or the other you have your data 150 print "initial guess for big-C (=10^{- little-c} ) ?" 160 input C 170 Orig_c=C 180 'note: C=0.1 thru 1.4 OK for cooked data (C=1); beyond that, underflow 190 Oldc=-1000 200 'loop here when starting with a new c 210 if abs(C-Oldc)<0.001*C then 650 220 if abs(C-Orig_c)>(Orig_c+1) then 670 230 N1=0:N2=0:N3=0:N4=0:N5=0:N6=0:N7=0:N8=0:N9=0:N10=0:N11=0:N12=0:N13=0:N14=0:N15=0:N16=0:N17=0:N18=0:N19=0 240 for I=1 to 8 250 W(I)=1/(1+C*X(I)) 260 N1=N1+1 270 N2=N2+W(I) 280 N3=N3+W(I)^2 290 N4=N4+Y1(I):N5=N5+Y2(I):N6=N6+Y3(I) 300 N7=N7+Y1(I)*W(I):N8=N8+Y2(I)*W(I):N9=N9+Y3(I)*W(I) 310 N10=N10+X(I)*W(I)^2 320 N11=N11+X(I)*W(I)^3 330 N12=N12+X(I)*Y1(I)*W(I)^2:N13=N13+X(I)*Y2(I)*W(I)^2:N14=N14+X(I)*Y3(I)*W(I)^2 340 N15=N15+X(I)^2*W(I)^3 350 N16=N16+X(I)^2*W(I)^4 360 N17=N17+X(I)^2*Y1(I)*W(I)^3:N18=N18+X(I)^2*Y2(I)*W(I)^3:N19=N19+X(I)^2*Y3(I)*W(I)^3 370 next I 380 'for I=1 to 8:print W(I),:next I 390 print 400 'print N1,N2,N3,N4,N5,N6 410 Det=N1*N3-N2*N2:if Det<0.0001 then print "warning: small denominator" 420 B1=(N3*N4-N2*N7)/Det 430 A1=(-N2*N4+N1*N7)/Det 440 B2=(N3*N5-N2*N8)/Det 450 A2=(-N2*N5+N1*N8)/Det 460 B3=(N3*N6-N2*N9)/Det 470 A3=(-N2*N6+N1*N9)/Det 480 Top=N10*A1:Bottom=N10*B1+2*N11*A1-N12 490 B1p=(N3*Top-N2*Bottom)/Det:' The "p" stands for "-prime" (derivative) 500 A1p=(-N2*Top+N1*Bottom)/Det 510 Top=N10*A2:Bottom=N10*B2+2*N11*A2-N13 520 B2p=(N3*Top-N2*Bottom)/Det 530 A2p=(-N2*Top+N1*Bottom)/Det 540 Top=N10*A3:Bottom=N10*B3+2*N11*A3-N14 550 B3p=(N3*Top-N2*Bottom)/Det 560 A3p=(-N2*Top+N1*Bottom)/Det 570 F=-A1*N12-A2*N13-A3*N14+(A1*B1+A2*B2+A3*B3)*N10+(A1^2+A2^2+A3^2)*N11 580 print "latest value of F(C) is ",F:'hoping for F(C)=0 590 print "latest value of C is ",C 600 Fp=2*(A1*N17+A2*N18+A3*N19)-(A1p*N12+A2p*N13+A3p*N14)-2*(A1*B1+A2*B2+A3*B3)*N15+(A1*B1p+A1p*B1+A2*B2p+A2p*B2+A3*B3p+A3p*B3)*N10-3*(A1^2+A2^2+A3^2)*N16+2*(A1*A1p+A2*A2p+A3*A3p)*N11 610 Oldc=C 620 C=C-F/Fp:'Newton's method: C(n+1)=C(n)-F(C(n))/F'(C(n)) 630 'end of loop 640 goto 200 650 print "well, that should be enough; I quit." 660 goto 680 670 print "this doesn't seem to be working; I give up" 680 'get here when tired of looping 690 print "Chosen values of parameters are C=";C 700 print "A1, B1, A2, B2, A3, B3=",A1,B1,A2,B2,A3,B3 710 end:' a variant of your sample dataset follows 720 data 0.205,8.51,1.98,2.36 730 data 0.223,8.51,1.96,2.34 740 data 0.253,8.54,1.96,2.31 750 data 0.283,8.64,1.96,2.28 760 data 0.323,9.21,2.03,2.24 770 data 0.394,9.21,2.03,2.24 780 data 0.491,9.54,2.06,2.18 790 data 0.591,9.63,2.06,2.12 ============================================================================== Date: Tue, 7 Feb 1995 09:46:28 +0000 (GMT) From: [Permission pending] To: Dave Rusin Subject: Re: program to find your C Excellent - that's wonderful, thanks. I'll give it a try as soon as I can. The data definately is of the form (eqn that I sent), because I can fit each individual dataset to that with R^2 > 0.95. [sig deleted] ============================================================================== Date: Tue, 14 Feb 95 00:38:52 CST From: rusin (Dave Rusin) To: [Permission pending] Subject: Re: program to find your C Did you ever succeed in fitting the data to curves of the sort you wanted? I know when I sent you the basic program it worked just fine for artificial data which was made to be a close fit, but it didn't succeed for the one dataset you posted. Did your other data fit the curves sufficiently well that you could adapt the program for your use? If not, what did you do instead? dave (always curious) ============================================================================== Date: Tue, 14 Feb 1995 11:07:04 +0000 (GMT) From: [Permission pending] To: Dave Rusin Subject: Re: program to find your C Hi, > Did you ever succeed in fitting the data to curves of the sort you wanted? I did (eventually), thanks. > I know when I sent you the basic program it worked just fine for artificial > data which was made to be a close fit, but it didn't succeed for the one > dataset you posted. Did your other data fit the curves sufficiently well > that you could adapt the program for your use? If not, what did you do > instead? > I couldn't get your program to work (thanks for doing it, btw). The rest of my data was actually worse than the set I posted (that's the problem with experimental data sometimes, I'm afraid). I managed to do it with Mathematica (thanks to someone's kind help), and that worked well. I tihnk it would have failed awfully with more than about 10 points (it was a crude but effective way), but for my data it worked fine. All I've got to do now is finish writing my thesis ... > dave (always curious) > Thanks again for your help [sig deleted - djr]