### Dawson - Matlab How To

 For the assignment of calculating various probability distributions with given forecast and observed data, I wrote a script with Matlab to handle the process of making counts from the data and calculating probabilities based on those counts. This initially required the forecast and observed data to be added as arrays in Matlab using the “Import Data” tool. After these data were added, the steps listed below were taken.The blocks of code show examples of how this work could be done with Matlab. These examples assume one set of forecasts and corresponding observations for a dichotomous (yes or no) event. Therefore, a 2x2 contingency table is used to calculate the respective probability distributions.1. Make a scatter plot of the joint distribution of forecasts and observationsThe figure command simply creates a new figure, rather than overwriting an existing oneUse the scatter function to make the scatter plotThe '.' defines the type of marker on the scatter plotUse the title, xlabel, and ylabel functions to define the title and axes labels for the plot`figure``scatter(forecasts(:),observations(:),'.');``title('Joint Distribution of Forecasts and Observations')``xlabel('Forecasts')``ylabel('Observations')`2. Create 2-D array to represent the contingency table (R-forecasts, C-observations)Use 'for loop' to parse through rows of forecast and observation data and establish counts for each element of the contingency tableUse ‘if statements’ to determine the row for the forecast and the column for the observationWithin the if statements, define row (i) and column (j) index for the appropriate forecast and observation categoriesAfter each iteration of the for loop, add 1 to the value in the appropriate contingency table cell`contingency_table(1:2,1:2) = 0;``for k=1:length(forecasts);``    if forecasts(k) == 1``        i = 1;``    else``        i = 2;``    end``    if observations(k) == 1``        j = 1;``    else``        j = 2;``    end``    contingency_table(i,j) = contingency_table(i,j) + 1;``end`3. Calculate joint distributionCreate joint distribution arrayShould have same dimensions as contingency tableDivide all elements in the contingency table by NN = number of forecast-observation pairs`joint(1:2,1:2) = 0;``joint(:,:) = contingency_table(:,:)/N;`4. Calculate marginal distributionsCreate forecast and observation count arraysShould have dimensions of number of rows for forecast and columns for observationsUse for loop to sum up the rows and columns of contingency table to make count arraysCreate marginal distribution arraysShould have same dimensions as count arraysUse for loop to parse through count arrays and divide each count by the total number of forecasts/observations`fcount(1:2,1) = 0;``xcount(1,1:2) = 0;``for i=1:length(contingency_table)``    fcount(i,:) = sum(contingency_table(i,:));``    xcount(:,i) = sum(contingency_table(:,i));``end`pf(1:2,1) = 0;px(1,1:2) = 0;for i=1:length(pf);    pf(i,1) = fcount(i,1)/length(forecasts);    px(1,i) = xcount(1,i)/length(observations);end5. Calculate conditional distributionsCreate conditional distribution arraysShould have same dimensions as contingency table/joint distributionUse nested for loops to count through marginal distributions of forecasts and observationsUse if statement to account for the possibility of dividing by a joint probability of zeroWhen p(x) and/or p(f) is not equal to zero, divide each element of the joint distribution by the appropriate marginal distribution`pfx(1:2,1:2) = 0;``pxf(1:2,1:2) = 0;``for i=1:length(pf);``    for j=1:length(px);``        if px(j) == 0``            pfx(i,j) = NaN;``        elseif pf(i) == 0``            pxf(i,j) = NaN;``        else``            pfx(i,j) = joint(i,j)/px(j);``            pxf(i,j) = joint(i,j)/pf(i);``        end``    end``end`6. Calculate forecast biasCalculate average value of all forecasts using mean functionCalculate average value of all observations using mean functionSubtract average observed value from average forecast value`avg_f = mean(forecasts(:));``avg_o = mean(observations(:));``bias = avg_f - avg_o;`