One of the most useful bits of number-crunching you can do with data is to calculate the probability distribution, in the earnest hope that it will be a reasonable fit for one of the recognised distributions such as the normal distribution. In this project I will write a small C library to calculate the normal distribution for a given data set.
Most people will have an intuitive understanding of the normal distribution even if they have no formal education in statistics. As an example let's think about people's heights: a few people are very short, rather more are fairly short, a lot are around average height, and then the numbers drop back down as we start to consider heights above average, eventually reaching zero.
If we draw a graph of heights and their respective frequencies we might expect to see something like this. (Note this is just fictitious data I have created for this project.)
The bars show actual data, whereas the red circles show the theoretical normal distribution for this particular data. As you can see there is a very good but not perfect fit, for example the frequency of 66 inches is 21, but the theoretical frequency assuming the data fits a normal distribution is 20.1516.
So how do we calculate, for example, 20.1516 from our data? Short answer: using the formula below, borrowed from Wikipedia. The variables used in this formula are:
The formula looks slightly intimidating but later on I'll implement it in C to calculate the normal values for each item of data. In fact that is at the core of this whole project - to iterate a set of data and apply the above formula to calculate the values plotted as red circles.
The Plan
Before starting to code let's consider our requirements and how we are going to implement them. Basically what we want to end up with is a data structure with a row for each different value in the original data, and columns for the following:
The value (eg. height)
The actual frequency of that value
The actual probability of that value
The theoretical frequency according to the normal distribution
The theoretical probability according to the normal distribution
For simplicity I am restricting the data to integers, each of which will be counted separately when calculating frequencies. A refinement would be to allow for real values (float or double in C) and to group these into class intervals for the purposes of calculating frequencies. A further refinement could be to extend the code to calculate other types of probability distribution.
Coding
Create a new folder and within it create the following empty files. You can download the source code as a zip or pull it from Github if you prefer.
data.h
data.c
normaldistribution.h
normaldistribution.c
main.c
ZIP File GitHub
In the real world we would be getting our data from a database or some kind of file but for the purposes of this project I have written a function to populate an array with hard-coded sample data. In a couple of earlier projects on this site I have used randomly generated sample data but here we need data which roughly follows a normal distribution so have resorted to hard coding it. I have on my list of future projects a program to generate random data following certain probability distributions but for now this will have to suffice.
In data.h enter the function prototype.
data.h
bool populatedata(int data[8]);
Now open data.c and enter the function body.
data.c
#include<stdio.h>
The normaldistribution.h header file contains two struct declarations. The central one is prob_dist which will hold an array of individual items, one for each discrete value, as well as overall totals. The other struct is prob_dist_item to hold values, frequencies and probabilities. We'll see exactly how this works later.
Also in the header file are the four function prototypes.
normaldistribution.h
#include<stdlib.h>
That's the preliminaries and infrastructure done so we can get stuck in to the central pieces of code. Open normaldistribution.c and enter or paste the following. I'll run through it in detail afterwards.
normaldistribution.c
#include<stdio.h>
Firstly we need to #define a couple of constants if they don't already exist - these are Pi (π) and e, or Euler's Number.
Next there are a couple of function prototypes for static functions which I'll describe when we get to them.
The normal_distribution_create function uses malloc to get enough memory for a prob_dist struct and then initializes the members to 0 or NULL.
Now we move on to normal_distribution_calculate, the central function of the whole project.
Firstly we'll declare a few variables for later use, and set pd->total_frequency to size.
Next we iterate the data. On each iteration if we haven't encountered the value before we add a new prob_dist_item to the array and set the frequency to 1. If we have already seen the value before it will already have a prob_dist_item with its associated frequency, so we just increment the frequency. To check which is the case I have used a separate static function called indexof. Along the way we also increment total and sumofsquares for later use.
We now need to sort the items. This is not strictly essential but users aren't going to want their numbers all jumbled up. I have used the qsort function from stdlib.h. I wrote a post about qsort a while ago so won't go into details here, but briefly it takes an array to sort, the size of the array, the size of an individual item, and a pointer to a comparator function which we must provide ourselves. How you implement this is up to you but it takes void pointers to the two items to compare and must return <1, 0 or >1 depending on whether the first item is less than, equal to or more than the second. Our version is compare_prob_dist_item.
We can now calculate the arithmetic mean, variance and standard deviation, before incrementing the recently created array of prob_dist_item structs to fill in the missing values.
For each discrete value we calculate the following three values:
probability - simply the frequency divided by size normal_probability - this uses the equation shown above but repeated here normal_frequency - this is the normal probability multiplied by the number of data items
We also increment the three totals, total_probability, total_normal_probability and total_normal_frequency.
Now let's look at the normal_distribution_print function. This prints out a heading and then iterates the prob_dist_item array, printing out the contents in a table.
The normal_distribution_free function is very simple; it just calls free on the prob_dist struct's array of prob_dist_item structs, and then on the prob_dist struct itself.
The main function creates an int array and then calls the various functions we have written. To recap and provide an overview of the entire process I'll list the steps:
Create an int array
Populate it with populatedata
Create a prob_dist struct
Calculate actual and normal distribution frequencies and probabilities with normal_distribution_calculate
Print these with normal_distribution_print
Free up dynamic memory with normal_distribution_free
main.c
#include<stdio.h>
We have now finished the coding so to compile and run the program run the following in your terminal.
Compile and Run
gcc main.c data.c normaldistribution.c -std=c11 -lm -o main
This is the program output. The normal probability and frequency totals don't quite add up to 1 and 138 respectively due to roundings.
Program Output
-----------------------
