/* nag_cp_stat (g02ecc) Example Program. * * Copyright 2002 Numerical Algorithms Group. * * Mark 7, 2002. * Mark 7b revised, 2004. */ #include #include #include #include #include #include int main(void) { /* Scalars */ double sigsq, tss; Integer exit_status, num_models, i, ii, j, m, n, nmod, pdx; NagError fail; Nag_OrderType order; /* Arrays */ char **model = 0; double *cp = 0, *rsq = 0, *rss = 0, *wtptr = 0, *x = 0, *y = 0; Integer *sx = 0, *mrank = 0, *nterms = 0; const char *var_names[] = { "DAY", "BOD", "TKN", "TS", "TVS", "COD" }; /* For iteration over model */ Integer model_index = 0; #ifdef NAG_COLUMN_MAJOR #define X(I, J) x[(J-1)*pdx + I - 1] order = Nag_ColMajor; #else #define X(I, J) x[(I-1)*pdx + J - 1] order = Nag_RowMajor; #endif INIT_FAIL(fail); exit_status = 0; printf("nag_cp_stat (g02ecc) Example Program Results\n"); /* Skip heading in data file */ scanf("%*[^\n] "); scanf("%ld%ld%*[^\n] ", &n, &m); /* Allocate memory */ if (!(x = NAG_ALLOC(n * m, double)) || !(y = NAG_ALLOC(n, double)) || !(sx = NAG_ALLOC(m, Integer))) { printf("Allocation failure\n"); exit_status = -1; goto END; } #ifdef NAG_COLUMN_MAJOR pdx = n; order = Nag_ColMajor; #else pdx = m; order = Nag_RowMajor; #endif for (i = 1; i <= n; ++i) { for (j = 1; j <= m; ++j) scanf("%lf", &X(i, j)); scanf("%lf%*[^\n] ", &y[i - 1]); } num_models = 1; for (j = 1; j <= m; ++j) { scanf("%ld", &sx[j - 1]); if (sx[j - 1] == 1) num_models <<= 1; } scanf("%*[^\n] "); /* Allocate memory */ if (!(model = NAG_ALLOC(num_models*m, char *)) || !(cp = NAG_ALLOC(num_models, double)) || !(rsq = NAG_ALLOC(num_models, double)) || !(rss = NAG_ALLOC(num_models, double)) || !(mrank = NAG_ALLOC(num_models, Integer)) || !(nterms = NAG_ALLOC(num_models, Integer))) { printf("Allocation failure\n"); exit_status = -1; goto END; } /* Calculate residual sums of squares using nag_all_regsn (g02eac) */ /* nag_all_regsn (g02eac). * Computes residual sums of squares for all possible linear * regressions for a set of independent variables */ nag_all_regsn(order, Nag_MeanInclude, n, m, x, pdx, var_names, sx, y, wtptr, &nmod, (const char**)model, rss, nterms, mrank, &fail); if (fail.code != NE_NOERROR) { printf("Error from nag_all_regsn (g02eac).\n%s\n", fail.message); exit_status = 1; goto END; } tss = rss[0]; sigsq = rss[nmod - 1] / (n - nterms[nmod - 1] - 1); /* nag_cp_stat (g02ecc). * Calculates R^2 and C_P values from residual sums of * squares */ nag_cp_stat(Nag_MeanInclude, n, sigsq, tss, nmod, nterms, rss, rsq, cp, &fail); if (fail.code != NE_NOERROR) { printf("Error from nag_cp_stat (g02ecc).\n%s\n", fail.message); exit_status = 1; goto END; } printf("\n"); printf("Number of CP RSQ MODEL\n"); printf("parameters\n"); for (i = 1; i <= nmod; ++i) { ii = nterms[i - 1]; printf(" %7ld%11.2f%8.4f ", ii, cp[i - 1], rsq[i - 1]); for (j = 1; j <= ii; ++j) printf("%-3.3s %s", model[model_index++], j%5 == 0 || j == 5?"\n":" "); printf("\n"); } END: if (model) NAG_FREE(model); if (cp) NAG_FREE(cp); if (rsq) NAG_FREE(rsq); if (rss) NAG_FREE(rss); if (x) NAG_FREE(x); if (y) NAG_FREE(y); if (sx) NAG_FREE(sx); if (mrank) NAG_FREE(mrank); if (nterms) NAG_FREE(nterms); return exit_status; }