Statistics::Smoothing::SGT - A Simple Good-Turing (SGT) smoothing implementation

  use Statistics::Smoothing::SGT
  my $sgt = new Statistics::Smoothing::SGT($frequencyClasses, $total);
  $sgt->calculateValues();
  $probabilities = $sgt->getProbabilities();
  $newFrequencies = $sgt->getNewFrequencies();
  $nBar = $sgt->getNBar();

    The Free Software Foundation, Inc.,
    59 Temple Place - Suite 330,
    Boston, MA  02111-1307, USA.

package Statistics::Smoothing::SGT;

# use strict, as we do not our variables to go haywire
use strict;

# use these for debugging purposes
use warnings;
use Data::Dumper;

our ($VERSION);

$VERSION = '2.1.2';

# constructor
sub new {
	my ($class, $frequencyClasses, $ngramTotal) = @_;
	my $rowsFound = scalar(keys(%{$frequencyClasses}));
	unless ($rowsFound >=5) {
	    die("At least 5 m/V(m) pairs must be provided for SGT.\nThe hash contains only $rowsFound rows.\n");
	}
	my $self = {
		frequencyClasses => $frequencyClasses,
		ngramTotal => $ngramTotal
	};
	bless($self, $class);
	return $self;
}

# method for calculating all Z(m)
sub calculateZ() {
	my ($self) = @_;
	my @zValues;
	my @equivalenceClasses;
	my @cardinalities;
	
	# iterate over frequencies for conversion into ascendingly ordered list
	foreach my $frequency (sort {$a <=> $b} keys(%{$self->{frequencyClasses}})) {
		# push frequency to array of equivalence classes
		push(@equivalenceClasses, $frequency);

		# push cardinality of this frequency class to array of cardinalities at corresponding
		# position
		push(@cardinalities, $self->{frequencyClasses}->{$frequency});
	}
	
	# save arrays for processing in other functions
	$self->{equivalenceClasses} = \@equivalenceClasses;
	$self->{cardinalities} = \@cardinalities;
	
	# calculate z value for m = 1
	push(@zValues, 2 * $self->{cardinalities}->[0] / ($self->{equivalenceClasses}->[1]));
	
	# iterate over equivalence classes and their respective cardinalities
	# for calculating z values for all m > 1 and m < max
	for (my $i = 1; $i < @{$self->{equivalenceClasses}} - 1; $i++) {
		push(@zValues, 2 * $self->{cardinalities}->[$i] / ($self->{equivalenceClasses}->[$i + 1] - $self->{equivalenceClasses}->[$i - 1]));
	}

	# calculate z value for m = max
	push(@zValues, 2 * $self->{cardinalities}->[@{$self->{cardinalities}} - 1] / ($self->{equivalenceClasses}->[@{$self->{equivalenceClasses}} - 1] - $self->{equivalenceClasses}->[@{$self->{equivalenceClasses}} - 2]));

	# return
	return \@zValues;
}

# method for getting first gap between frequencies
sub getGap {
	my ($self) = @_;
	my $gapAt;
	my $buffer = 0;

	# iterate over equivalence classes
	for (my $i = 0; $i <= @{$self->{equivalenceClasses}}; $i++) {
		# check if gap between current value and previous
		# one is bigger than 1
		if (!defined $self->{equivalenceClasses}->[$i] || $buffer + 1 < $self->{equivalenceClasses}->[$i]) {
			# if so, we have found our gap and can skip the
			# remaining iterations
			$gapAt = $buffer;
			last;
		}

		# write current value to buffer for checking next value
		$buffer = $self->{equivalenceClasses}->[$i];
	}
	
	# return gap
	return $gapAt;
}

# method for deducing slope and intersection of a logarithmised
# exponential function with its best fitting linear function
# (i.e. linear regression)
sub logBestFit { 
	my ($self, $Xaxis, $Yaxis) = @_;
	my $XYs = 0;
	my $Xsquares = 0;
	my $meanX = 0;
	my $meanY = 0;
	my $rows = scalar @{$Xaxis};

	# calculate log mean values for x and y
	for (my $i = 0; $i < $rows; $i++) {
		$meanX += log($Xaxis->[$i]);
		$meanY += log($Yaxis->[$i]);
	}
	$meanX /= $rows;
	$meanY /= $rows;

	# calculate slope and intersection
	for (my $i = 0; $i < $rows; $i++) {
		$XYs += ((log($Xaxis->[$i]) - $meanX) * (log($Yaxis->[$i]) - $meanY));
		$Xsquares += (log($Xaxis->[$i]) - $meanX) ** 2;
	}
	my $slope = $XYs / $Xsquares;
	my $intersection = $meanY - $slope * $meanX;
	
	# return slope and intersection
	return ($slope, $intersection);
}

# public method for calculating SGT-smoothed values
sub calculateValues() {
	my ($self) = @_;
	my $gapAt;
	my $slope;
	my $intersection;
	my %xValues;
	my %yValues;
	my $zValues;
	my %mFound;
	my %mStar;
	my $newFrequencies;
	my $m;
	my $xFlag;
	my $diff;
	my $i;

	# get z values
	$zValues = $self->calculateZ();

	# get first frequency gap
	$gapAt = $self->getGap();

	# perform linear regression
	($slope, $intersection) = $self->logBestFit($self->{equivalenceClasses}, $zValues);
	
	# calculcate y function values
    for ($i = 0; $i <= $self->{equivalenceClasses}->[@{$self->{equivalenceClasses}} - 1]; $i++) {
		$yValues{$i} = ($i + 1) * exp($slope * (log($i + 2) - log($i + 1)));
    }

	# calculate x value only, if last element has not yet
	# been reached, x values for large classes do not matter
	# anyway
	for($i = 0; $i < @{$self->{equivalenceClasses}} - 2; $i++) {
		$xValues{$self->{equivalenceClasses}->[$i]} = $self->{equivalenceClasses}->[$i + 1] * $self->{cardinalities}->[$i + 1] / $self->{cardinalities}->[$i];
	}

	# build hash for all existing m (=equivalence class) and V(m) (= cardinality)
	for ($i = 0; $i < @{$self->{equivalenceClasses}}; $i++){
		$mFound{$self->{equivalenceClasses}->[$i]} = $self->{cardinalities}->[$i]; 
	}

	# iterate over equivalence classes in order to determine
	# which function to use
	for ($i = 1; $i <= $self->{equivalenceClasses}->[@{$self->{equivalenceClasses}} - 1]; $i++) {
		# check, if m (=equivalence class) and V(m) (= cardinality) exist for this index
		unless ($mFound{$i}) {
			next;
		}
		
		# if y function has not been used yet and first gap between frequencies
		# has not yet been reached
		if (!$xFlag && $i < $gapAt) {
			# calculate distance value between x and y function for which y
			# function is to be used
			$diff = (sqrt (($self->{equivalenceClasses}->[$i + 1] ** 2) *
					($self->{cardinalities}->[$i + 1] / $self->{cardinalities}->[$i] ** 2) *
					(1 + $self->{cardinalities}->[$i + 1] / $self->{cardinalities}->[$i]))) * 1.65;

			# if difference between x and y value is smaller than the difference
			# value calculated above
			if (abs($xValues{$i} - $yValues{$i}) < $diff) {
				# use y function from now on
				$mStar{$i} = $yValues{$i};
				$xFlag = 1;
			} else {
				# use x function
				$mStar{$i} = $xValues{$i};
			}
		} else {
			# use y function from now on
			$mStar{$i} = $yValues{$i};
			$xFlag = 1;
		}
		
		# calculate new total (i.e. nBar);
		$self->{nBar} += $mStar{$i} * $mFound{$i};
	}

	# save new frequencies
	$self->{newFrequencies} = \%mStar;

	# calculate the probability for unseen events (i.e. pZero)
	$self->{probabilities}->{0} = $self->{cardinalities}->[0] / $self->{ngramTotal};

	# calculate all other probabilities
	foreach $m (keys(%mStar)) {
		$self->{probabilities}->{$m} = (1 - $self->{probabilities}->{0}) * ($mStar{$m} / $self->{nBar});
	}
}

# public method for getting new frequencies
sub getNewFrequencies() {
	my ($self) = @_;

	# return equivalence classes
	return $self->{newFrequencies};
}

# public method for getting equivalence classes
sub getEquivalenceClasses() {
	my ($self) = @_;

	# return equivalence classes
	return $self->{equivalenceClasses};
}

# public method for getting probabilities
sub getProbabilities() {
	my ($self) = @_;

	# return probabilities
	return $self->{probabilities};
}

# public method for getting nBar
sub getNBar() {
	my ($self) = @_;

	# return nBar
	return $self->{nBar};
}

1;

__END__