当前位置： 文江博客话题详情

使用 nls 再现 PROC NLIN 输出

发布于 2024-11-23 22:15:23 字数 10586 浏览 1 评论 0原文

我想在 R 中实现以下 SAS 代码。

$$N_e = N_o{1-exp[\frac{(d+bN_o)(T_h N_e - T)}{(1+c N_o)}]}$$

其中 $b>0$，$c\geq 0$ ，$T_h>0$，并且$T=72$。

SAS代码是

 DATA NOTONECT;
  INPUT N0  REP FATE    NE  N02 N03 PROPEAT;
  DATALINES;
5   1   0   0   25  125 0
5   2   0   1   25  125 0.2
5   3   0   1   25  125 0.2
5   4   0   2   25  125 0.4
5   5   0   2   25  125 0.4
5   6   0   2   25  125 0.4
5   7   0   2   25  125 0.4
5   8   0   3   25  125 0.6
7   1   0   0   49  343 0
7   2   0   0   49  343 0
7   3   0   1   49  343 0.14286
7   4   0   1   49  343 0.14286
7   5   0   2   49  343 0.28571
7   6   0   2   49  343 0.28571
7   7   0   2   49  343 0.28571
7   8   0   3   49  343 0.42857
10  1   0   1   100 1000    0.1
10  2   0   1   100 1000    0.1
10  3   0   2   100 1000    0.2
10  4   0   2   100 1000    0.2
10  5   0   3   100 1000    0.3
10  6   0   3   100 1000    0.3
10  7   0   3   100 1000    0.3
10  8   0   4   100 1000    0.4
10  9   0   7   100 1000    0.7
15  1   0   1   225 3375    0.06667
15  2   0   1   225 3375    0.06667
15  3   0   3   225 3375    0.2
15  4   0   3   225 3375    0.2
15  5   0   4   225 3375    0.26667
15  6   0   5   225 3375    0.33333
15  7   0   5   225 3375    0.33333
15  8   0   5   225 3375    0.33333
20  1   0   3   400 8000    0.15
20  2   0   4   400 8000    0.2
20  3   0   7   400 8000    0.35
20  4   0   7   400 8000    0.35
20  5   0   8   400 8000    0.4
20  6   0   8   400 8000    0.4
20  7   0   9   400 8000    0.45
20  8   0   11  400 8000    0.55
25  1   0   4   625 15625   0.16
25  2   0   5   625 15625   0.2
25  3   0   6   625 15625   0.24
25  4   0   7   625 15625   0.28
25  5   0   9   625 15625   0.36
25  6   0   9   625 15625   0.36
25  7   0   13  625 15625   0.52
25  8   0   14  625 15625   0.56
30  1   0   5   900 27000   0.16667
30  2   0   8   900 27000   0.26667
30  3   0   10  900 27000   0.33333
30  4   0   11  900 27000   0.36667
30  5   0   11  900 27000   0.36667
30  6   0   12  900 27000   0.4
30  7   0   14  900 27000   0.46667
30  8   0   20  900 27000   0.66667
45  1   0   4   2025    91125   0.08889
45  2   0   7   2025    91125   0.15556
45  3   0   8   2025    91125   0.17778
45  4   0   10  2025    91125   0.22222
45  5   0   11  2025    91125   0.24444
45  6   0   14  2025    91125   0.31111
45  7   0   15  2025    91125   0.33333
45  8   0   19  2025    91125   0.42222
60  1   0   9   3600    216000  0.15
60  2   0   14  3600    216000  0.23333
60  3   0   14  3600    216000  0.23333
60  4   0   16  3600    216000  0.26667
60  5   0   18  3600    216000  0.3
60  6   0   21  3600    216000  0.35
60  7   0   24  3600    216000  0.4
60  8   0   26  3600    216000  0.43333
80  1   0   7   6400    512000  0.0875
80  2   0   11  6400    512000  0.1375
80  3   0   12  6400    512000  0.15
80  4   0   15  6400    512000  0.1875
80  5   0   17  6400    512000  0.2125
80  6   0   12  6400    512000  0.15
80  7   0   21  6400    512000  0.2625
80  8   0   23  6400    512000  0.2875
100 1   0   7   10000   1000000 0.07
100 2   0   8   10000   1000000 0.08
100 3   0   10  10000   1000000 0.1
100 4   0   11  10000   1000000 0.11
100 5   0   15  10000   1000000 0.15
100 6   0   24  10000   1000000 0.24
100 7   0   26  10000   1000000 0.26
100 8   0   33  10000   1000000 0.33
;

PROC NLIN DATA=NOTONEC2;
 PARMS BHAT= 0.001 0.01 0.1 /* initial parameter estimates */
 CHAT= 0.001 0.01 0.1 
 DHAT= 0 THHAT=3.0; 
 BOUNDS BHAT>0,CHAT>=0,
 THHAT>0; /* parameter bounds */
 T=72; /* experimental period in H */
 X=NE; /* initial predicted value */ 
 A=(DHAT+BHAT*N0)/(1+CHAT*N0); /* expression for A */
 /* define the implicit function */
 C1=EXP(-A*T); /* components of the implicit function */
 C2=A*THHAT;
 H=N0*C1*EXP(C2*X)+X-N0; /* the implicit function */
 ITER=0; /* iterations for Newton's method */
 /* Newton's method employed to find predicted number eaten */
 DO WHILE(ABS(H)>0.0001 AND ITER<50); /* stop criteria for Newton's method */
      X=X-H/(N0*C1*C2*EXP(C2*X)+1); /* new predicted value */
      H=N0*C1*EXP(C2*X)+X-N0; /* new value of implicit function */
      ITER=ITER+1; /* iteration counter */
 END;
MODEL NE=X; /* model for nonlinear least squares */
Run;

SAS输出是

                                       The NLIN Procedure
                                  Dependent Variable NE

                                       Grid Search
                                                                  Sum of
                    BHAT        CHAT        DHAT       THHAT     Squares

                 0.00100     0.00100           0      3.0000      3714.9
                  0.0100     0.00100           0      3.0000     47602.6
                  0.1000     0.00100           0      3.0000      100896
                 0.00100      0.0100           0      3.0000      3118.3
                  0.0100      0.0100           0      3.0000     29395.6
                  0.1000      0.0100           0      3.0000     97841.3
                 0.00100      0.1000           0      3.0000      1804.1
                  0.0100      0.1000           0      3.0000      6625.5
                  0.1000      0.1000           0      3.0000     53940.3



                                   The NLIN Procedure
                                  Dependent Variable NE
                                  Method: Gauss-Newton

                                     Iterative Phase
                                                                       Sum of
             Iter        BHAT        CHAT        DHAT       THHAT     Squares

                0     0.00100      0.1000           0      3.0000      1804.1
                1    0.000396      0.0176     0.00190      3.2240      1774.5
                2    0.000327     0.00578     0.00206      3.5044      1762.4
                3    0.000307     0.00193     0.00208      3.6507      1754.8
                4    0.000302    0.000793     0.00208      3.7006      1752.3
                5    0.000299    0.000265     0.00208      3.7246      1751.1
                6    0.000299    0.000123     0.00208      3.7312      1750.7
                7    0.000298    0.000039     0.00208      3.7351      1750.6
                8    0.000298    0.000013     0.00208      3.7363      1750.5
                9    0.000298    5.224E-6     0.00208      3.7367      1750.5
               10    0.000298    2.547E-7     0.00208      3.7369      1750.5
               11    0.000298    6.092E-8     0.00208      3.7369      1750.5
               12    0.000298    3.59E-10     0.00208      3.7369      1750.5
               13    0.000298    6.37E-11     0.00208      3.7369      1750.5
               14    0.000298    -287E-13     0.00208      3.7369      1750.5
               15    0.000298    -864E-13     0.00208      3.7369      1750.5
               16    0.000298    -955E-13     0.00208      3.7369      1750.5
               17    0.000298    -983E-13     0.00208      3.7369      1750.5
               18    0.000298    -992E-13     0.00208      3.7369      1750.5
               19    0.000298    -997E-13     0.00208      3.7369      1750.5
               20    0.000298    -999E-13     0.00208      3.7369      1750.5
               21    0.000298      -1E-10     0.00208      3.7369      1750.5
               22    0.000298      -1E-10     0.00208      3.7369      1750.5
               23    0.000298      -1E-10     0.00208      3.7369      1750.5
               24    0.000298      -1E-10     0.00208      3.7369      1750.5
               25    0.000298      -1E-10     0.00208      3.7369      1750.5
               26    0.000298      -1E-10     0.00208      3.7369      1750.5
               27    0.000298      -1E-10     0.00208      3.7369      1750.5
               28    0.000407           0    0.000884      4.0885      1733.8
               29    0.000416           0    0.000825      4.0645      1733.2
               30    0.000415           0    0.000832      4.0640      1733.2
               31    0.000416           0    0.000832      4.0640      1733.2


            NOTE: Convergence criterion met.


                                   Estimation Summary

                          Method                   Gauss-Newton
                          Iterations                         31
                          Subiterations                     108
                          Average Subiterations        3.483871
                          R                             1.51E-6

                                   The NLIN Procedure

                                   Estimation Summary

                          PPC(DHAT)                    0.000057
                          RPC(DHAT)                    0.000663
                          Object                       2.84E-10
                          Objective                    1733.203
                          Observations Read                  89
                          Observations Used                  89
                          Observations Missing                0


                  NOTE: An intercept was not specified for this model.

                                          Sum of        Mean               Approx
        Source                    DF     Squares      Square    F Value    Pr > F

        Model                      3      9526.8      3175.6     157.57    <.0001
        Error                     86      1733.2     20.1535
        Uncorrected Total         89     11260.0


                                   Approx       Approximate 95%
     Parameter      Estimate    Std Error      Confidence Limits     Label

     BHAT           0.000416     0.000228    -0.00004    0.000869
     CHAT                  0            0           0           0
     DHAT           0.000832      0.00382    -0.00676     0.00842
     THHAT            4.0640       0.3575      3.3534      4.7747
     Bound0            107.9        315.6      -509.1       725.0    0 <= CHAT


                             Approximate Correlation Matrix
                           BHAT            CHAT            DHAT           THHAT

          BHAT        1.0000000               .      -0.8928900       0.7190236
          CHAT         .                      .        .               .
          DHAT       -0.8928900               .       1.0000000      -0.5474038
          THHAT       0.7190236               .      -0.5474038       1.0000000

我不知道如何在nls中使用牛顿法。如果有人帮助我解决这个问题，我将不胜感激。提前致谢。

原文

I'd like to implement the following SAS code in R.

$$N_e = N_o{1-exp[\frac{(d+bN_o)(T_h N_e - T)}{(1+c N_o)}]}$$

where $b>0$, $c\geq 0$, $T_h>0$, and $T=72$.

The SAS code is

 DATA NOTONECT;
  INPUT N0  REP FATE    NE  N02 N03 PROPEAT;
  DATALINES;
5   1   0   0   25  125 0
5   2   0   1   25  125 0.2
5   3   0   1   25  125 0.2
5   4   0   2   25  125 0.4
5   5   0   2   25  125 0.4
5   6   0   2   25  125 0.4
5   7   0   2   25  125 0.4
5   8   0   3   25  125 0.6
7   1   0   0   49  343 0
7   2   0   0   49  343 0
7   3   0   1   49  343 0.14286
7   4   0   1   49  343 0.14286
7   5   0   2   49  343 0.28571
7   6   0   2   49  343 0.28571
7   7   0   2   49  343 0.28571
7   8   0   3   49  343 0.42857
10  1   0   1   100 1000    0.1
10  2   0   1   100 1000    0.1
10  3   0   2   100 1000    0.2
10  4   0   2   100 1000    0.2
10  5   0   3   100 1000    0.3
10  6   0   3   100 1000    0.3
10  7   0   3   100 1000    0.3
10  8   0   4   100 1000    0.4
10  9   0   7   100 1000    0.7
15  1   0   1   225 3375    0.06667
15  2   0   1   225 3375    0.06667
15  3   0   3   225 3375    0.2
15  4   0   3   225 3375    0.2
15  5   0   4   225 3375    0.26667
15  6   0   5   225 3375    0.33333
15  7   0   5   225 3375    0.33333
15  8   0   5   225 3375    0.33333
20  1   0   3   400 8000    0.15
20  2   0   4   400 8000    0.2
20  3   0   7   400 8000    0.35
20  4   0   7   400 8000    0.35
20  5   0   8   400 8000    0.4
20  6   0   8   400 8000    0.4
20  7   0   9   400 8000    0.45
20  8   0   11  400 8000    0.55
25  1   0   4   625 15625   0.16
25  2   0   5   625 15625   0.2
25  3   0   6   625 15625   0.24
25  4   0   7   625 15625   0.28
25  5   0   9   625 15625   0.36
25  6   0   9   625 15625   0.36
25  7   0   13  625 15625   0.52
25  8   0   14  625 15625   0.56
30  1   0   5   900 27000   0.16667
30  2   0   8   900 27000   0.26667
30  3   0   10  900 27000   0.33333
30  4   0   11  900 27000   0.36667
30  5   0   11  900 27000   0.36667
30  6   0   12  900 27000   0.4
30  7   0   14  900 27000   0.46667
30  8   0   20  900 27000   0.66667
45  1   0   4   2025    91125   0.08889
45  2   0   7   2025    91125   0.15556
45  3   0   8   2025    91125   0.17778
45  4   0   10  2025    91125   0.22222
45  5   0   11  2025    91125   0.24444
45  6   0   14  2025    91125   0.31111
45  7   0   15  2025    91125   0.33333
45  8   0   19  2025    91125   0.42222
60  1   0   9   3600    216000  0.15
60  2   0   14  3600    216000  0.23333
60  3   0   14  3600    216000  0.23333
60  4   0   16  3600    216000  0.26667
60  5   0   18  3600    216000  0.3
60  6   0   21  3600    216000  0.35
60  7   0   24  3600    216000  0.4
60  8   0   26  3600    216000  0.43333
80  1   0   7   6400    512000  0.0875
80  2   0   11  6400    512000  0.1375
80  3   0   12  6400    512000  0.15
80  4   0   15  6400    512000  0.1875
80  5   0   17  6400    512000  0.2125
80  6   0   12  6400    512000  0.15
80  7   0   21  6400    512000  0.2625
80  8   0   23  6400    512000  0.2875
100 1   0   7   10000   1000000 0.07
100 2   0   8   10000   1000000 0.08
100 3   0   10  10000   1000000 0.1
100 4   0   11  10000   1000000 0.11
100 5   0   15  10000   1000000 0.15
100 6   0   24  10000   1000000 0.24
100 7   0   26  10000   1000000 0.26
100 8   0   33  10000   1000000 0.33
;

PROC NLIN DATA=NOTONEC2;
 PARMS BHAT= 0.001 0.01 0.1 /* initial parameter estimates */
 CHAT= 0.001 0.01 0.1 
 DHAT= 0 THHAT=3.0; 
 BOUNDS BHAT>0,CHAT>=0,
 THHAT>0; /* parameter bounds */
 T=72; /* experimental period in H */
 X=NE; /* initial predicted value */ 
 A=(DHAT+BHAT*N0)/(1+CHAT*N0); /* expression for A */
 /* define the implicit function */
 C1=EXP(-A*T); /* components of the implicit function */
 C2=A*THHAT;
 H=N0*C1*EXP(C2*X)+X-N0; /* the implicit function */
 ITER=0; /* iterations for Newton's method */
 /* Newton's method employed to find predicted number eaten */
 DO WHILE(ABS(H)>0.0001 AND ITER<50); /* stop criteria for Newton's method */
      X=X-H/(N0*C1*C2*EXP(C2*X)+1); /* new predicted value */
      H=N0*C1*EXP(C2*X)+X-N0; /* new value of implicit function */
      ITER=ITER+1; /* iteration counter */
 END;
MODEL NE=X; /* model for nonlinear least squares */
Run;

The SAS output is

                                       The NLIN Procedure
                                  Dependent Variable NE

                                       Grid Search
                                                                  Sum of
                    BHAT        CHAT        DHAT       THHAT     Squares

                 0.00100     0.00100           0      3.0000      3714.9
                  0.0100     0.00100           0      3.0000     47602.6
                  0.1000     0.00100           0      3.0000      100896
                 0.00100      0.0100           0      3.0000      3118.3
                  0.0100      0.0100           0      3.0000     29395.6
                  0.1000      0.0100           0      3.0000     97841.3
                 0.00100      0.1000           0      3.0000      1804.1
                  0.0100      0.1000           0      3.0000      6625.5
                  0.1000      0.1000           0      3.0000     53940.3



                                   The NLIN Procedure
                                  Dependent Variable NE
                                  Method: Gauss-Newton

                                     Iterative Phase
                                                                       Sum of
             Iter        BHAT        CHAT        DHAT       THHAT     Squares

                0     0.00100      0.1000           0      3.0000      1804.1
                1    0.000396      0.0176     0.00190      3.2240      1774.5
                2    0.000327     0.00578     0.00206      3.5044      1762.4
                3    0.000307     0.00193     0.00208      3.6507      1754.8
                4    0.000302    0.000793     0.00208      3.7006      1752.3
                5    0.000299    0.000265     0.00208      3.7246      1751.1
                6    0.000299    0.000123     0.00208      3.7312      1750.7
                7    0.000298    0.000039     0.00208      3.7351      1750.6
                8    0.000298    0.000013     0.00208      3.7363      1750.5
                9    0.000298    5.224E-6     0.00208      3.7367      1750.5
               10    0.000298    2.547E-7     0.00208      3.7369      1750.5
               11    0.000298    6.092E-8     0.00208      3.7369      1750.5
               12    0.000298    3.59E-10     0.00208      3.7369      1750.5
               13    0.000298    6.37E-11     0.00208      3.7369      1750.5
               14    0.000298    -287E-13     0.00208      3.7369      1750.5
               15    0.000298    -864E-13     0.00208      3.7369      1750.5
               16    0.000298    -955E-13     0.00208      3.7369      1750.5
               17    0.000298    -983E-13     0.00208      3.7369      1750.5
               18    0.000298    -992E-13     0.00208      3.7369      1750.5
               19    0.000298    -997E-13     0.00208      3.7369      1750.5
               20    0.000298    -999E-13     0.00208      3.7369      1750.5
               21    0.000298      -1E-10     0.00208      3.7369      1750.5
               22    0.000298      -1E-10     0.00208      3.7369      1750.5
               23    0.000298      -1E-10     0.00208      3.7369      1750.5
               24    0.000298      -1E-10     0.00208      3.7369      1750.5
               25    0.000298      -1E-10     0.00208      3.7369      1750.5
               26    0.000298      -1E-10     0.00208      3.7369      1750.5
               27    0.000298      -1E-10     0.00208      3.7369      1750.5
               28    0.000407           0    0.000884      4.0885      1733.8
               29    0.000416           0    0.000825      4.0645      1733.2
               30    0.000415           0    0.000832      4.0640      1733.2
               31    0.000416           0    0.000832      4.0640      1733.2


            NOTE: Convergence criterion met.


                                   Estimation Summary

                          Method                   Gauss-Newton
                          Iterations                         31
                          Subiterations                     108
                          Average Subiterations        3.483871
                          R                             1.51E-6

                                   The NLIN Procedure

                                   Estimation Summary

                          PPC(DHAT)                    0.000057
                          RPC(DHAT)                    0.000663
                          Object                       2.84E-10
                          Objective                    1733.203
                          Observations Read                  89
                          Observations Used                  89
                          Observations Missing                0


                  NOTE: An intercept was not specified for this model.

                                          Sum of        Mean               Approx
        Source                    DF     Squares      Square    F Value    Pr > F

        Model                      3      9526.8      3175.6     157.57    <.0001
        Error                     86      1733.2     20.1535
        Uncorrected Total         89     11260.0


                                   Approx       Approximate 95%
     Parameter      Estimate    Std Error      Confidence Limits     Label

     BHAT           0.000416     0.000228    -0.00004    0.000869
     CHAT                  0            0           0           0
     DHAT           0.000832      0.00382    -0.00676     0.00842
     THHAT            4.0640       0.3575      3.3534      4.7747
     Bound0            107.9        315.6      -509.1       725.0    0 <= CHAT


                             Approximate Correlation Matrix
                           BHAT            CHAT            DHAT           THHAT

          BHAT        1.0000000               .      -0.8928900       0.7190236
          CHAT         .                      .        .               .
          DHAT       -0.8928900               .       1.0000000      -0.5474038
          THHAT       0.7190236               .      -0.5474038       1.0000000

I don't know how to use Newton method in nls. I'd highly appreciate if someone help me to figure out this problem. Thanks in advance.

如果你对这篇内容有疑问，欢迎到本站社区发帖提问 参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论， 你可以免费 注册 一个本站的账号。

评论（1）

面犯桃花 2024-11-30 22:15:23

我没有和你完全一样的方法，但我想答案应该是可以的。

设置计算

s0 <- list(b=0.1,c=0.1,d=0,Th=3)
X <- read.table("rogersdat.txt",header=TRUE)

预期食用数量的起始条件函数：

predfun <- function(b,c,d,Th,Te,N0,debug=FALSE) {
  a <- (d+b*N0)/(1+c*N0)
  r <- N0 - (1/(a*Th))*lambertW(a*Th*N0*exp(a*(Th*N0-Te)))
  if (debug) cat(mean(a),b,c,d,Th,mean(r),"\n")
  r
}

绘制数据（只是为了确保）：

with(X,plot(Ne~N0))
## check starting value
lines(1:100,with(c(s0,X),predfun(b,c,d,Th,Te=72,N0=1:100)))

拟合：

library(emdbook)
n1 <- nls(Ne~predfun(b,c,d,Th,Te=72,N0),data=X,
    lower=c(1e-6,1e-6,-Inf,1e-6),algorithm="port",
    start=list(b=0.1,c=0.1,d=0.002,Th=3))
summary(n1)

Formula: Ne ~ predfun(b, c, d, Th, Te = 72, N0)

Parameters:
    Estimate Std. Error t value Pr(>|t|)   
b  0.0004155  0.0008745   0.475  0.63591   
c  0.0000010  0.0657237   0.000  0.99999   
d  0.0008318  0.0067374   0.123  0.90203   
Th 4.0639937  1.4665102   2.771  0.00686 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 4.516 on 85 degrees of freedom

Algorithm "port", convergence message: relative convergence (4) 

confint.default(n1)

          2.5 %      97.5 %
b  -0.001298523 0.002129547
c  -0.128815083 0.128817083
d  -0.012373153 0.014036799
Th  1.189686504 6.938300956

估计值匹配得很好，置信区间则不然（无论如何，这些类型的置信区间在边界上有点冒险...... ）

我实际上建议使用二项式误差的最大似然估计会更好一些。

pdat <- data.frame(N0=1:100,Ne=predict(n1,newdata=data.frame(N0=1:100)))
with(pdat,lines(N0,Ne,col=2))
library(ggplot2)
ggplot(X,aes(x=N0,y=Ne))+stat_sum(aes(size=factor(..n..)),alpha=0.5)+theme_bw()+
  geom_line(data=pdat,colour="red")

拟合值图

I haven't done it exactly the same way as you, but I think the answer should be OK.

Set up starting conditions

s0 <- list(b=0.1,c=0.1,d=0,Th=3)
X <- read.table("rogersdat.txt",header=TRUE)

Function for computed expected number eaten:

predfun <- function(b,c,d,Th,Te,N0,debug=FALSE) {
  a <- (d+b*N0)/(1+c*N0)
  r <- N0 - (1/(a*Th))*lambertW(a*Th*N0*exp(a*(Th*N0-Te)))
  if (debug) cat(mean(a),b,c,d,Th,mean(r),"\n")
  r
}

Plot the data (just to make sure):

with(X,plot(Ne~N0))
## check starting value
lines(1:100,with(c(s0,X),predfun(b,c,d,Th,Te=72,N0=1:100)))

Fit:

library(emdbook)
n1 <- nls(Ne~predfun(b,c,d,Th,Te=72,N0),data=X,
    lower=c(1e-6,1e-6,-Inf,1e-6),algorithm="port",
    start=list(b=0.1,c=0.1,d=0.002,Th=3))
summary(n1)

Formula: Ne ~ predfun(b, c, d, Th, Te = 72, N0)

Parameters:
    Estimate Std. Error t value Pr(>|t|)   
b  0.0004155  0.0008745   0.475  0.63591   
c  0.0000010  0.0657237   0.000  0.99999   
d  0.0008318  0.0067374   0.123  0.90203   
Th 4.0639937  1.4665102   2.771  0.00686 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 4.516 on 85 degrees of freedom

Algorithm "port", convergence message: relative convergence (4) 

confint.default(n1)

          2.5 %      97.5 %
b  -0.001298523 0.002129547
c  -0.128815083 0.128817083
d  -0.012373153 0.014036799
Th  1.189686504 6.938300956

Estimates match pretty well, confidence intervals don't (these types of confidence intervals are a little dicey on the boundary anyway ...)

I would actually suggest that maximum likelihood estimation with binomial errors would be a little better.

pdat <- data.frame(N0=1:100,Ne=predict(n1,newdata=data.frame(N0=1:100)))
with(pdat,lines(N0,Ne,col=2))
library(ggplot2)
ggplot(X,aes(x=N0,y=Ne))+stat_sum(aes(size=factor(..n..)),alpha=0.5)+theme_bw()+
  geom_line(data=pdat,colour="red")

plot of fitted values

回复 原文
~没有更多了~

关于作者

落在眉间の轻吻

暂无简介

0 文章
0 评论
22 人气
发私信

相关话题

更多

热门标签

操作系统 程序设计 IT运维 Linux系统管理 JavaScript 服务器应用 solaris C/C++ PHP Shell BSD Vue.js aix Oracle Python HTML 系统管理 HTML5 CSS 前端
更多

推荐作者

離殇

文章 0 评论 0

小姐丶请自重

文章 0 评论 0

Aik

文章 0 评论 0

国产ˉ祖宗

文章 0 评论 0

猥琐帝

文章 0 评论 0

半仙

文章 0 评论 0

更多

友情链接

文江博客
    我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站，即表示您同意使用 Cookies 和您的相关数据。
    原文