📄 ocn_g12.h
字号:
/*------------------------------------------------------------------------------*
* File Name: OCN_g12.h *
* Creation: SDB 5/22/2001 *
* Purpose: Origin C Header for NAG functions *
* Copyright (c) OriginLab Corp. 2001 *
* All Rights Reserved *
* *
* Modification Log: *
*------------------------------------------------------------------------------*/
#ifndef _O_NAG_G12_H
#define _O_NAG_G12_H
//#importdll "ONAG" // NAG DLL prepared by OriginLab
#pragma dll(ONAG)
#include <NAG\nag_types.h>
/* begin proto */
/** g12aac
computes the Kaplan-Meier,(or product-limit),estimates of
survival probabilities for a sample of failure times.
Example 1:
Assume we have "Data1" Worksheet which has 4 columns and 10 rows.
The first column contain 10 data, and we will put the result,
the time in the second column, survival probability in the third
column, and standard deviation in fourth column. We also have
"Data2" Worksheet which has 2 columns and 10 rows, these two columns
contain 10 data in each column.
int n = 10, nd;
//Attach six Datasets to these 6 columns
Dataset dt("data1",0);
Dataset dtp("data1",1);
Dataset dp("data1",2);
Dataset dpsig("data1",3);
Dataset difreq("data2",0);
Dataset dic("data2",1);
dt.SetSize(10);
dtp.SetSize(10);
dp.SetSize(10);
dpsig.SetSize(10);
difreq.SetSize(10);
dic.SetSize(10);
//Dataset cannot be the parameter of function, but vector can be.
vector t = dt, tp = dtp, p = dp, psig = dpsig;
vector<int> ifreq = difreq, ic = dic;
nag_prod_limit_surviv_fn(n, t, ic, ifreq, &nd, tp, p, psig);
//Put the result to worksheet "Data1" the second, third and fourth column.
dtp = tp;
dp = p;
dpsig = psig;
Example 2:
The remission times for a set of 21 leukemia patients at 18 distinct
time points are read in and the Kaplan-Meier estimate computed and
printed. For further details see Gross and Clark (1975), page 242.
void test_nag_prod_limit_surviv_fn()
{
double psig[18];
double p[18];
double t[18] = {6.0, 6.0, 7.0, 9.0, 10.0, 10.0, 11.0, 13.0, 16.0, 17.0,
19.0, 20.0, 22.0, 23.0, 25.0, 32.0, 34.0, 35.0};
double tp[18];
int i, n = 18, nd;
int ifreq[18] = {1, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1};
int ic[18] = {1, 0, 0, 1, 0, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1};
nag_prod_limit_surviv_fn(n, t, ic, ifreq, &nd, tp, p, psig);
printf("nd = %d \n",nd);
printf("\n Time Survival Standard\n");
printf(" probability deviation\n\n");
for (i = 0; i < nd; ++i)
printf(" %6.1f%10.3f %10.3f\n", tp[i], p[i], psig[i]);
}
The output is following:
nd = 7
Time Survival Standard
probability deviation
6.0 0.857 0.076
7.0 0.807 0.087
10.0 0.753 0.096
13.0 0.690 0.107
16.0 0.627 0.114
22.0 0.538 0.128
23.0 0.448 0.135
Return:
This function returns NAG error code, 0 if no error.
11: On entry, n must not be less than 2: n=_value_.
604: On entry, ic[_value_]=_value_. The censor code for an observation must be either 0 or 1.
605: On entry, freq[_value_]=_value_. The value of frequency for an observation must be = 0.
73: Memory allocation failed.
74: An internal error has occurred in this function. Check the function call and array sizes. If the function call is correct, please consult NAG for assistance.
successfully call of the nag_prod_limit_surviv_fn function.
*/
int nag_prod_limit_surviv_fn(
int n, //the number of failure and censored times give in t.
const double t[], //the failure and censored times;these need not be ordered.
const int ic[], //contains the censoring code of the ith observation.
const int freq[], //indicates whether frequencies are provided for each failure and censored time point.
int *nd, //the number of distinct failure times.
double tp[], //contains the ith ordered distinct failure time.
double p[], //contains the Kaplan-Meier estimate of the survival probability.
double psig[] //contains an estimate of the standard deviation.
); // General block design or completely randomized design.
/** g12bac
returns parameter estimates and other statistics that are associated with the Cox proportional
hazards model for fixed covariates
Example1:
Assume we have "Data1" Worksheet which is double type and has 4 columns
and 42 rows. The first two columns contain 42 data in each column, and we
will put the result, the time in the third column, survival function in
the fourth column. We have "Data2" Worksheet which is integer type and has
2 columns and 42 rows. The first column contains 42 data.
double b[2], cov[3], dev, df;
double omega[42], res[42], sc[2], se[2];
double tol = 5e-5, v[336], y[42], offset[42];
int i, i_1, ip = 1, ip1, iprint = 0, irank, isi[42], sz[1] = {1};
int j, m = 1, maxit = 20, n = 42, nd, ndmax = 42, ns = 0, tdsur = 1, tdv = 8;
double c_b42 = 0;
int c_0 = 0;
string filename = "c:\\test.txt";
ip1 = ip +1;
i_1 = n;
//Attach five Datasets to these five columns
Dataset dt("data1",0);
Dataset dz("data1",1);
Dataset dtp("data1",2);
Dataset dsur("data1",3);
Dataset dic("data2",0);
dt.SetSize(42);
dz.SetSize(42);
dtp.SetSize(42);
dsur.SetSize(42);
//Dataset cannot be the parameter of function, but vector can be.
vector t = dt, z = dz, tp = dtp, sur = dsur;
vector<int> ic = dic;
for(i =0; i < n; i++)
{
y[i] = 1.0 -1.0*ic[i];
v[i+250] = log(t[i]);
offset[i] = log(t[i]);
}
nag_glm_poisson(Nag_Log, Nag_MeanInclude, n, z, m, m, sz, ip1, y, NULL, offset,
c_b42, &dev, &df, b, &irank, se, cov, v, tdv, tol, maxit, c_0,
filename, c_b42);
for(i = 0; i < ip; i++)
b[i] = b[i+1];
if(irank != ip +1)
printf("WARNING: covariate not of full rank\n");
nag_surviv_cox_model(n, m, ns, z, m, sz, ip, t, ic, NULL, isi, &dev, b,
se, sc, cov, res, &nd, tp, sur, tdsur, ndmax, tol,
maxit, iprint, filename);
//put the result to Worksheet "data1", the third time is time and fourth column
//is suvival function.
dtp = tp;
dsur = sur;
Example2:
The data are the remission times for two groups of leukemia patients(see Gross
and Clark(1975)p242). A dummy variable indicates which group they come from. An
initial estimate is computed using the exponential model and then the Cox proportional
hazard model is fitted and parameter estimates and the survival function are printed
void test_nag_surviv_cox_model()
{
double b[2], cov[3], dev, df;
double omega[42], res[42], sc[2], se[2];
double sur[42];
double t[42] = {1, 1, 2, 2, 3, 4, 4, 5, 5, 8, 8, 8, 8, 11, 11, 12, 12,
15, 17, 22, 23, 6, 6, 6, 7, 10, 13, 16, 22, 23, 6, 9, 10, 11, 17, 19,
20, 25, 32, 32, 34, 35};
double z[42] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
double tol, tp[42], v[336], y[42], offset[42];
int i, i_1, ip, ip1, iprint = 0, irank, isi[42], sz[1];
int ic[42] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
int j, m = 1, maxit = 20, n = 42, nd, ndmax = 42, ns = 0, tdsur = 1, tdv;
double c_b42 = 0;
int c_0 = 0;
//Get the origin folder path.
char sztemp[100];
LT_get_str("%Y", sztemp, 100);
//Define a filename and working path.
string str1, filename;
str1 = sztemp;
filename = str1 + "text.txt";
printf("%s\n", filename);
sz[0] =1;
ip =1;
ip1 = ip +1;
tol = 5e-5;
i_1 = n;
tdv = 8;
for(i =0; i < n; i++)
{
y[i] = 1.0 -1.0*ic[i];
v[i+250] = log(t[i]);
offset[i] = log(t[i]);
}
nag_glm_poisson(Nag_Log, Nag_MeanInclude, n, z, m, m, sz, ip1, y, NULL, offset,
c_b42, &dev, &df, b, &irank, se, cov, v, tdv, tol, maxit, c_0,
filename, c_b42);
for(i = 0; i < ip; i++)
b[i] = b[i+1];
if(irank != ip +1)
printf("WARNING: covariate not of full rank\n");
nag_surviv_cox_model(n, m, ns, z, m, sz, ip, t, ic, NULL, isi, &dev, b,
se, sc, cov, res, &nd, tp, sur, tdsur, ndmax, tol,
maxit, iprint, filename);
printf("\n");
printf("Parameter Estimate");
printf(" Standard Error");
printf("\n");
for(i = 0; i < ip; i++)
printf("%6ld %8.4f %8.4f\n",i + 1, b[i], se[i]);
printf("\n");
printf("Deviance = %13.4f\n",dev);
printf("\n");
printf(" Time Survivor Function");
printf("\n");
ns = 1;
for(i = 0; i < nd; i++)
{
printf("%10.0f",tp[i]);
for(j = 0; j < ns; j++)
{
printf(" %8.4f",sur[i * tdsur +j]);
if( (j + 1) % 3 == 0)
printf("\n");
}
printf("\n");
}
}
The output is following:
C:\sftest\Gr394\text.txt
Parameter Estimate Standard Error
1 -1.5091 0.4096
Deviance = 172.7592
Time Survivor Function
1 0.9640
2 0.9264
3 0.9065
4 0.8661
5 0.8235
6 0.7566
7 0.7343
8 0.6506
10 0.6241
11 0.5724
12 0.5135
13 0.4784
15 0.4447
16 0.4078
17 0.3727
22 0.2859
23 0.1908
Return:
This function returns NAG error code, 0 if no error.
11: On entry, m must not be less than 1: m = <value>. On entry, n must not be less than 2: n = <value>. On entry, ns must not be less than 0: ns = <value>. On entry, tdsur must not be less than 1: tdsur = <value>. On entry, max_iter must not be less than 0: max_iter = <value>.
17: On entry, tdsur = <value> while ns = <value>. These parameters must satisfy tdsur >= ns. On entry, tdz = <value> while m = <value>. These parameters must satisfy tdz >= m.
117: On entry, tol = <value>, machine precision (X02AJC) = <value>. Constraint: tol >= 10 * machine precision.
78: Cannot open file outfile for appending.
79: Cannot close file outfile.
103: On entry, sz[<value>] =<value>. Constrant: sz[<value>] >= 0. On entry, ic[<value>] = <value>. Constraint: ic[<value>] = 0 or 1. On entry, isi[<value>] = <value>. Constraint: 0 <= isi[<value>] <= ns.
607: On entry, ip = <value> and the number of non zero values of sz = <value>. Constraint: ip = the number of non-zero values of sz.
608: On entry, the number of values of sz[i] > 0 is <value>, n = <value> and excluded observations with isi[i] = 0 is <value>. Constraint: the number of values of non-zero sz must be less than n - excluded observations.
609: On entry, ndmax is = <value> while the output value of nd = <value>. Constraint: ndmax >= nd.
610: The matrix of second partial derivatives is singular. Try different starting values or include fewer covariates.
611: Overflow has been detected. Try different starting values.
612: Convergence has not been achieved in max_iter iterations. The progress towards convergence can be examined by using by setting iprint to >= 1. Any non _convergence decreasing the deviance maybe due to a linear combination of covariates being monotonic with time. Full results are returned.
613: In the current iteration 10 step halvings have been performed withour decreasing the deviance from the previous iteration. Convergence is assumed.
73: Memory allocation failed.
74: An internal error has occurred in this function. Check the function call and array sizes. If the function call is correct, please consult NAG for assistance.
successfully call of the nag_surviv_cox_model function.
*/
int nag_surviv_cox_model(
int n, //the number of data points
int m, //the number of covariates
int ns, //the number of strata
const double z[], //the covariates
int tdz, //the second dimension of array z;
const int sz[], //indicates which subset of covariates is to be included
int ip, //the number of covariates included as indicated by sz
const double t[], //the vector of n failure censoring times
const int ic[], //the status of the individual at time t given in t
const double omega[], //if an offset is set, it contain value \omega[i], otherwise it is NULL;
const int isi[], //the stratum indicator
double *dev, //the deviance
double b[], //initial estimates of the covariate coefficient parameters
double se[], //the asymptotic standard error of the estimate
double sc[], //the value of the score function
double cov[], //the variance-covariance matirx of the parameter estimates
double res[], //the residuals
int *nd, //the number of distinct failure times
double tp[], //the distinct failure time
double sur[], //the estimated survival function
int tdsur, //the second dimension of the array
int ndmax, //the first dimension of the array
double tol, //the accuracy required for the estimation
int max_iter, //the maximun number of iterations
int iprint, //indicates if the printing information on the iteration is required
const char *outfile //the name of the file into which information is to be output.
);
#endif //!_O_NAG_G12_H⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -