/*-------------------------------------------------------------------- File : spcdiv.c Contents: compute specificity divergence for simple possibility distribution example Author : Christian Borgelt History : 20.10.1999 file created 06.01.2002 weighted sum of possibilities added 08.01.2002 estimation from data table added --------------------------------------------------------------------*/ #include #include #include /*-------------------------------------------------------------------- Preprocessor Definitions --------------------------------------------------------------------*/ #define LN_2 0.69314718055994530942 /* ln(2) */ #define log2(x) (log(x)/LN_2) #define TPLCNT (sizeof(tab)/sizeof(TUPLE)) #define max(x,y) (((x) > (y)) ? (x) : (y)) #define min(x,y) (((x) < (y)) ? (x) : (y)) /*-------------------------------------------------------------------- Type Definitions --------------------------------------------------------------------*/ typedef struct { /* --- data tuple --- */ int a, b, c; /* attribute values */ int cnt; /* number of occurrences */ } TUPLE; /* (tuple) */ /*-------------------------------------------------------------------- Constants --------------------------------------------------------------------*/ const TUPLE tab[] = { /* --- data tuples */ { 0, 0, 0, 60 }, { 0, 0, 1, 30 }, { 0, 0, 2, 10 }, { 0, 1, 0, 20 }, { 0, 1, 1, 10 }, { 0, 1, 2, 10 }, { 0, 2, 0, 10 }, { 0, 2, 1, 20 }, { 0, 2, 2, 30 }, { 1, 0, 0, 70 }, { 1, 0, 1, 40 }, { 1, 0, 2, 10 }, { 1, 1, 0, 10 }, { 1, 1, 1, 10 }, { 1, 1, 2, 10 }, { 1, 2, 0, 10 }, { 1, 2, 1, 70 }, { 1, 2, 2, 50 }, { 2, 0, 0, 10 }, { 2, 0, 1, 10 }, { 2, 0, 2, 10 }, { 2, 1, 0, 30 }, { 2, 1, 1, 60 }, { 2, 1, 2, 10 }, { 2, 2, 0, 10 }, { 2, 2, 1, 10 }, { 2, 2, 2, 10 }, { 3, 0, 0, 10 }, { 3, 0, 1, 10 }, { 3, 0, 2, 10 }, { 3, 1, 0, 30 }, { 3, 1, 1, 50 }, { 3, 1, 2, 10 }, { 3, 2, 0, 10 }, { 3, 2, 1, 60 }, { 3, 2, 2, 50 }, { 0,-1, 1, 10 }, { 0, 0,-1, 20 }, { 0, 1,-1, 10 }, { 0, 2,-1, 10 }, { 1, 0,-1, 20 }, { 1, 2,-1, 10 }, { 2, 0,-1, 10 }, { 2, 1,-1, 10 }, { 3, 1,-1, 10 }, { 3, 2,-1, 10 } }; /* -1: unknown value */ /*-------------------------------------------------------------------- Global Variables --------------------------------------------------------------------*/ double pi_A[4]; /* marginal distrib. on attribute A */ double pi_B[3]; /* marginal distrib. on attribute B */ double pi_C[3]; /* marginal distrib. on attribute C */ double pi_AB[4][3]; /* marginal distrib. on AxB */ double pi_BC[3][3]; /* marginal distrib. on BxC */ double pi_AC[4][3]; /* marginal distrib. on AxC */ double pi_ABC[4][3][3]; /* full possibility distribution */ double pi[4][3][3]; /* graphical model distribution */ /*-------------------------------------------------------------------- Functions --------------------------------------------------------------------*/ double spcdiv (void) { /* --- specificity divergence */ int a, b, c; /* attribute values */ int n1, n2; /* size of alpha-cuts */ double alpha; /* alpha level (loop variable) */ double delta = 0.001; /* step width for alpha level */ double d = 0; /* specificity divergence */ for (alpha = delta; alpha <= 1; alpha += delta) { n1 = n2 = 0; for (a = 4; --a >= 0; ) { for (b = 3; --b >= 0; ) { for (c = 3; --c >= 0; ) { if (pi_ABC[a][b][c] >= alpha) n1++; if (pi [a][b][c] >= alpha) n2++; } /* compute the sizes of the */ } /* alpha-cuts of the two */ } /* possibility distributions */ if ((n1 > 0) && (n2 > 0)) d += delta *(log2(n2) -log2(n1)); } return d; } /* spcdiv() */ /*------------------------------------------------------------------*/ double wgtsum (void) { /* --- weighted sum of poss. degrees */ int i; /* loop variables */ int a, b, c; /* attribute values */ int ae, be, ce; /* ditto, for loops */ double s = 0, t; /* weighted sum of poss. degrees */ for (i = TPLCNT; --i >= 0; ) { t = 0; if (tab[i].a >= 0) { ae = tab[i].a; a = ae+1; } else { ae = 0; a = 4; } while (--a >= ae) { if (tab[i].b >= 0) { be = tab[i].b; b = be+1; } else { be = 0; b = 3; } while (--b >= be) { if (tab[i].c >= 0) { ce = tab[i].c; c = ce+1; } else { ce = 0; c = 3; } while (--c >= ce) if (pi[a][b][c] > t) t = pi[a][b][c]; } } s += tab[i].cnt *t; } return s; } /* wgtsum() */ /*------------------------------------------------------------------*/ int main (void) { int i; int a, b, c; int ae, be, ce; /* --- compute distribution --- */ for (i = TPLCNT; --i >= 0; ) { if (tab[i].a >= 0) { ae = tab[i].a; a = ae+1; } else { ae = 0; a = 4; } while (--a >= ae) { if (tab[i].b >= 0) { be = tab[i].b; b = be+1; } else { be = 0; b = 3; } while (--b >= be) { if (tab[i].c >= 0) { ce = tab[i].c; c = ce+1; } else { ce = 0; c = 3; } while (--c >= ce) pi_ABC[a][b][c] += tab[i].cnt; } } } for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi_ABC[a][b][c] /= 1000; /* --- compute marginal distributions --- */ for (a = 4; --a >= 0; ) { pi_A[a] = 0; for (b = 3; --b >= 0; ) { pi_AB[a][b] = 0; for (c = 3; --c >= 0; ) pi_AB[a][b] = max(pi_AB[a][b], pi_ABC[a][b][c]); pi_A[a] = max(pi_A[a], pi_AB[a][b]); } } for (b = 3; --b >= 0; ) { pi_B[b] = 0; for (c = 3; --c >= 0; ) { pi_BC[b][c] = 0; for (a = 4; --a >= 0; ) pi_BC[b][c] = max(pi_BC[b][c], pi_ABC[a][b][c]); pi_B[b] = max(pi_B[b], pi_BC[b][c]); } } for (c = 3; --c >= 0; ) { pi_C[c] = 0; for (a = 4; --a >= 0; ) { pi_AC[a][c] = 0; for (b = 3; --b >= 0; ) pi_AC[a][c] = max(pi_AC[a][c], pi_ABC[a][b][c]); pi_C[c] = max(pi_C[c], pi_AC[a][c]); } } /* --- print distributions --- */ printf("AxB|C=c1: AxB|C=c2: AxB|C=c3:\n"); for (b = 3; --b >= 0; ) { for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_ABC[a][b][0]); printf(" "); for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_ABC[a][b][1]); printf(" "); for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_ABC[a][b][2]); printf("\n"); } printf("\nAxB: BxC: AxC:\n"); for (b = 3; --b >= 0; ) { for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_AB[a][b]); printf(" "); for (c = 0; c < 3; c++) printf("%3.0f ", 1000 *pi_BC[b][c]); printf(" "); for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_AC[a][b]); printf("\n"); } printf("\nA: "); for (a = 0; a < 4; a++) printf("%3.0f ", 1000 *pi_A[a]); printf("\nB: "); for (b = 0; b < 3; b++) printf("%3.0f ", 1000 *pi_B[b]); printf("\nC: "); for (c = 0; c < 3; c++) printf("%3.0f ", 1000 *pi_C[c]); printf("\n\n"); /* --- specificity divergence --- */ printf("model spcdiv wgtsum\n"); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_A[a], min(pi_B[b], pi_C[c])); printf("A B C: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_AB[a][b], pi_C[c]); printf("A-B C: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_A[a], pi_BC[b][c]); printf("A B-C: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_AC[a][c], pi_B[b]); printf("A-C B: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_AB[a][b], pi_BC[b][c]); printf("A-B-C: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_AB[a][b], pi_AC[a][c]); printf("B-A-C: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = min(pi_AC[a][c], pi_BC[b][c]); printf("A-C-B: %.3f %4.1f\n", spcdiv(), wgtsum()); for (a = 4; --a >= 0; ) for (b = 3; --b >= 0; ) for (c = 3; --c >= 0; ) pi[a][b][c] = pi_ABC[a][b][c]; printf("A-B-C-A: %.3f %4.1f\n", spcdiv(), wgtsum()); return 0; } /* main() */