1 #include <stddef.h>
2 #include <string.h>
3 #include <stdio.h>
4 #include <math.h>
5 #include <Rmath.h>
6 #include <R.h>
7 #include "vector.h"
8 #include "subroutines.h"
9 #include "rand.h"
10 #include "bayes.h"
11 #include "sample.h"
12 
13 /* Prediction for Nonparametric Model for 2x2 Tables */
preDP(double * pdmu,double * pdSigma,int * pin_samp,int * pin_draw,int * pin_dim,int * verbose,double * pdStore)14 void preDP(
15 	   double *pdmu,
16 	   double *pdSigma,
17 	   int *pin_samp,
18 	   int *pin_draw,
19 	   int *pin_dim,
20 	   int *verbose,    /* 1 for output monitoring */
21 	   double *pdStore
22 	   ){
23 
24   /* some integers */
25   int n_samp = *pin_samp;    /* sample size */
26   int n_draw = *pin_draw;    /* sample size of survey data */
27   int n_dim = *pin_dim;      /* dimension */
28 
29   double *mu = doubleArray(n_dim);                /* The mean */
30   double *Wstar = doubleArray(n_dim);
31   double **Sigma = doubleMatrix(n_dim, n_dim);    /* The covariance matrix */
32 
33   /* misc variables */
34   int i, j, k, main_loop;   /* used for various loops */
35   int itemp = 0;
36   int itempM = 0;
37   int itempS = 0;
38   int progress = 1, itempP = ftrunc((double) n_draw/10);
39 
40   /* get random seed */
41   GetRNGstate();
42 
43   for(main_loop=0; main_loop<n_draw; main_loop++){
44     for(i=0; i<n_samp; i++) {
45       for (j=0;j<n_dim;j++) {
46 	mu[j] = pdmu[itempM++];
47 	for (k=j;k<n_dim;k++) {
48 	  Sigma[j][k] = pdSigma[itempS++];
49 	  Sigma[k][j] = Sigma[j][k];
50 	}
51       }
52       rMVN(Wstar, mu, Sigma, n_dim);
53       for (j=0; j<n_dim; j++)
54 	pdStore[itemp++] = exp(Wstar[j])/(1+exp(Wstar[j]));
55     }
56     if (*verbose)
57       if (itempP == main_loop) {
58         Rprintf("%3d percent done.\n", progress*10);
59         itempP+=ftrunc((double) n_draw/10); progress++;
60         R_FlushConsole();
61       }
62     R_CheckUserInterrupt();
63   }
64 
65   if(*verbose)
66     Rprintf("100 percent done.\n");
67 
68   /** write out the random seed **/
69   PutRNGstate();
70 
71   /* Freeing the memory */
72   free(mu);
73   free(Wstar);
74   FreeMatrix(Sigma,n_dim);
75 
76 } /* main */
77 
78