File: goodturing.c

package info (click to toggle)
bow 19991122-4
  • links: PTS
  • area: main
  • in suites: woody
  • size: 2,544 kB
  • ctags: 2,987
  • sloc: ansic: 38,660; lisp: 1,072; makefile: 594; perl: 492; yacc: 149; sh: 91
file content (228 lines) | stat: -rw-r--r-- 6,125 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
/* Simple Good-Turing smoothing */

/* Copyright (C) 1997, 1998, 1999 Andrew McCallum

   Written by:  Kamal Nigam <knigam@cs.cmu.edu> and 
                Rich Caruana <caruana@jprc.com>

   This file is part of the Bag-Of-Words Library, `libbow'.

   This library is free software; you can redistribute it and/or
   modify it under the terms of the GNU Library General Public License
   as published by the Free Software Foundation, version 2.
   
   This library is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
   Library General Public License for more details.

   You should have received a copy of the GNU Library General Public
   License along with this library; if not, write to the Free Software
   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111, USA */


/*
 *
 *
 *	Simple Good-Turing Frequency Estimator
 *
 *
 *	Geoffrey Sampson, with help from Miles Dennis
 *
 *	School of Cognitive and Computing Sciences
 *	University of Sussex
 *	England
 *
 *	27 June 1995
 *
 *
 *	Takes a set of (frequency, frequency-of-frequency) pairs, and 
 *	applies the Simple Good-Turing technique (devised by William A. Gale 
 *	of AT&T Bell Labs and described in Gale & Sampson, "Good-Turing
 *	Frequency Estimation Without Tears") for estimating 
 *	the probabilities corresponding to the observed frequencies, 
 *	and P.0, the joint probability of all unobserved species.
 *	Code is in ANSI standard C.
 *
 *	The program is written to take input from "stdin" and send output
 *	to "stdout"; redirection can be used to take input from and
 *	send output to permanent files.
 *
 *	The input file should be a series of lines separated by newline
 *	characters, where all nonblank lines contain two positive integers
 *	(an observed frequency, followed by the frequency of that frequency)
 *	separated by whitespace.  (Blank lines are ignored.)
 *	The lines should be in ascending order of frequency.
 *
 *	No checks are made for linearity; the program simply assumes that the 
 *	requirements for using the SGT estimator are met.
 *
 *	The output is a series of lines each containing an integer followed  
 *	by a probability (a real number between zero and one), separated by a 
 *	tab.  In the first line, the integer is 0 and the real number is the 
 *	estimate for P.0.  In subsequent lines, the integers are the  
 *	successive observed frequencies, and the reals are the estimated  
 *	probabilities corresponding to those frequencies.
 *
 *	No warranty is given as to absence of bugs.
 *
 *
 */
 



#include <bow/libbow.h>
#include <argp/argp.h>
#include <stdio.h>
#include <math.h>
#include <ctype.h>
#include <stdlib.h>
#include <string.h>

#define TRUE	1
#define FALSE	0
#define MAX_ROWS	100

static int r[MAX_ROWS], n[MAX_ROWS];
static double Z[MAX_ROWS], log_r[MAX_ROWS], log_Z[MAX_ROWS], 
  rStar[MAX_ROWS], p[MAX_ROWS];

static long int bigN;
static int rows;
static double PZero, bigNprime;
static double slope, intercept;

static double sq(double x)
{
  return(x * x);
}
 
static void findBestFit(void)
{
  double XYs, Xsquares, meanX, meanY;
  double sq(double);
  int i;
 	
  XYs = Xsquares = meanX = meanY = 0.0;
  for (i = 0; i < rows; ++i)
    {
      meanX += log_r[i];
      meanY += log_Z[i];
    }
  meanX /= rows;
  meanY /= rows;
  for (i = 0; i < rows; ++i)
    {
      XYs += (log_r[i] - meanX) * (log_Z[i] - meanY);
      Xsquares += sq(log_r[i] - meanX);
    }
  slope = XYs / Xsquares;
  intercept = meanY - slope * meanX;
}
 	
static double smoothed(int i)
{
  return(exp(intercept + slope * log(i)));
}
 	
static int row(int i)
{
  int j = 0;
 	
  while (r[j] < i && j < rows)
    ++j;
  return((j < rows && r[j] == i) ? j : -1);
}


/*
        - terms 0-(k-1) will be smoothed
        - freq must be of length k, and the last element must 
          be a legit count (SGT needs k count to smooth term (k-1))
        - simple_good_turing does not use freq[0]
        - disc must be of length k-1 or greater
        - arrays freq and disc have element 0 correspond to terms
          with 0 observed frequency, element 1 terms with frequency
          1, ...
        - because we are smoothing just terms 0-(k-1), probabilities 
          will not sum to one without external normalization because 
          terms >(k+1) where not smoothed
	  */
int simple_good_turing (int length, int *freq, double *disc)
{
  int i, j, next_n;
  double k, x, y;
  int indiffValsSeen = FALSE;
  int row(int);
  void findBestFit(void);
  double smoothed(int);
  double sq(double);
  void showEstimates(void);

  assert(length < MAX_ROWS);

  for (i=1; i<length; i++)
    {
      r[i-1] = i;
      n[i-1] = freq[i];
    }
  rows = length-1;

#if 0
  for (i=0; i<length; i++)
    printf ("%d %d\n", r[i], n[i]);
#endif

  bigN = 0;
  for (j = 0; j < rows; ++j)
    bigN += r[j] * n[j];
  PZero = n[row(1)] / (double) bigN;
  for (j = 0; j < rows; ++j)
    {
      i = (j == 0 ? 0 : r[j - 1]);
      if (j == rows - 1)
	k = (double) (2 * r[j] - i);
      else
	k = (double) r[j + 1];
      Z[j] = 2 * n[j] / (k - i);
      log_r[j] = log(r[j]);
      log_Z[j] = log(Z[j]);
    }
  findBestFit();
  for (j = 0; j < rows; ++j)
    {
      y = (r[j] + 1) * smoothed(r[j] + 1) / smoothed(r[j]);
      if (row(r[j] + 1) < 0)
	indiffValsSeen = TRUE;
      if (! indiffValsSeen)
	{
	  x = (r[j] + 1) * (next_n = n[row(r[j] + 1)]) / 
	    (double) n[j];
	  if (fabs(x - y) <= 1.96 * sqrt(sq(r[j] + 1.0) *
					 next_n / (sq((double) n[j])) * 
					 (1 + next_n / (double) n[j])))
	    indiffValsSeen = TRUE;
	  else
	    rStar[j] = x;
	}
      if (indiffValsSeen)
	rStar[j] = y;
    }
  bigNprime = 0.0;
  for (j = 0; j < rows; ++j)
    bigNprime += n[j] * rStar[j];
  for (j = 0; j < rows; ++j)
    p[j] = (1 - PZero) * rStar[j] / bigNprime;

  disc[0] = PZero;
  for (j=1; j<length; j++)
    {
      bow_verbosify(bow_progress, "(%d %f) ", j, p[j-1]);
      disc[j] = p[j-1] * (double) bigN / (double) j;
    }
  bow_verbosify(bow_progress, "\n");
  

  return(0);
}