<jupyter_start><jupyter_text>Exercise: Build a simple logistic regression modelIn this exercise, we'll fit a simple logistic regression model that will try to predict the chance of an avalanche. Recall that logistic regression fits an s-shaped curve to data, rather than a straight line, and we can use this to calculate a probability of a binary outcome. Data visualizationLet's start this exercise by loading in and having a look at our data:<jupyter_code>import pandas
!pip install statsmodels
!wget https://raw.githubusercontent.com/MicrosoftDocs/mslearn-introduction-to-machine-learning/main/graphing.py
!wget https://raw.githubusercontent.com/MicrosoftDocs/mslearn-introduction-to-machine-learning/main/Data/avalanche.csv

#Import the data from the .csv file
dataset = pandas.read_csv('avalanche.csv', delimiter="\t")

#Let's have a look at the data
dataset<jupyter_output>Requirement already satisfied: statsmodels in /anaconda/envs/py37_default/lib/python3.7/site-packages (0.12.2)
Requirement already satisfied: scipy>=1.1 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from statsmodels) (1.5.2)
Requirement already satisfied: numpy>=1.15 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from statsmodels) (1.19.2)
Requirement already satisfied: patsy>=0.5 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from statsmodels) (0.5.1)
Requirement already satisfied: pandas>=0.21 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from statsmodels) (1.2.0)
Requirement already satisfied: python-dateutil>=2.7.3 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from pandas>=0.21->statsmodels) (2.8.1)
Requirement already satisfied: pytz>=2017.3 in /anaconda/envs/py37_default/lib/python3.7/site-packages (from pandas>=0.21->statsmodels) (2020.5)
Requirement already satisfied: six in /anaconda/envs/py37_default/lib/pyth[...]<jupyter_text>Data ExplorationThe `avalanche` field is our target. A value of `1` means that an avalanche did occur at the conditions described by the features, whereas a value of `0` means no avalanche hapenned. Since our targets can only be `0` or `1` we call this a *binary classification* model.Now let's plot the relationships between each feature and the target values. That helps us understand which features are more likely to influence the results:<jupyter_code>import graphing # custom graphing code. See our GitHub repo for details

graphing.box_and_whisker(dataset, label_x="avalanche", label_y="surface_hoar", show=True)
graphing.box_and_whisker(dataset, label_x="avalanche", label_y="fresh_thickness", show=True)
graphing.box_and_whisker(dataset, label_x="avalanche", label_y="weak_layers", show=True)
graphing.box_and_whisker(dataset, label_x="avalanche", label_y="no_visitors")<jupyter_output><empty_output><jupyter_text>We can notice that:- For `fresh_thickness` the outcomes are very similar. This means that variations in their values aren't strongly correlated with the results.- Variations in values for `weak_layers` and `no_visitors`, seem to correlate with a larger number of `avalanche` results, and thus we should assign more importance to these features.The differences between avalanche and non-avalanche days are small and there isn't one clear driver of issues. Weak layers looks like a good starting point as it is related to the widest variation in results. Building a simple logistic regression modelWe will now built and train a model to predict the chance of an avalanche happening based __solely__ on the number of weak layers of snow:<jupyter_code># Here we import a function that splits datasets according to a given ratio
from sklearn.model_selection import train_test_split

# Split the dataset in an 70/30 train/test ratio. 
train, test = train_test_split(dataset, test_size=0.3, random_state=2)
print(train.shape)
print(test.shape)<jupyter_output>(766, 8)
(329, 8)<jupyter_text>OK, lets train our model, using the `train` dataset we've just created (notice that `weak_layers` will be the only feature used to determine the outcome):<jupyter_code>import statsmodels.formula.api as smf
import graphing # custom graphing code. See our GitHub repo for details

# Perform logistic regression.
model = smf.logit("avalanche ~ weak_layers", train).fit()

print("Model trained")<jupyter_output>Optimization terminated successfully.
         Current function value: 0.631451
         Iterations 5
Model trained<jupyter_text>After training, we can print a model summary with very detailed information:<jupyter_code>print(model.summary())<jupyter_output>Logit Regression Results                           
==============================================================================
Dep. Variable:              avalanche   No. Observations:                  766
Model:                          Logit   Df Residuals:                      764
Method:                           MLE   Df Model:                            1
Date:                Fri, 27 Aug 2021   Pseudo R-squ.:                 0.07898
Time:                        02:15:33   Log-Likelihood:                -483.69
converged:                       True   LL-Null:                       -525.17
Covariance Type:            nonrobust   LLR p-value:                 8.395e-20
===============================================================================
                  coef    std err          z      P>|z|      [0.025      0.975]
-------------------------------------------------------------------------------
Intercept      -0.8586      0.147     -5.856     [...]<jupyter_text>Notice that the positive coefficient for `weak_layers` means that a higher value means a higher likelihood for an avalanche.  Using our modelWe can now use our trained model to make predictions and estimate probabilities.Let's pick the first four occurrences in our `test` set and print the probability of an avalanche for each one of them:<jupyter_code># predict to get a probability

# get first 3 samples from dataset
samples = test["weak_layers"][:4]

# use the model to get predictions as possibilities
estimated_probabilities = model.predict(samples)

# Print results for each sample
for sample, pred in zip(samples,estimated_probabilities):
    print(f"A weak_layer with value {sample} yields a {pred * 100:.2f}% chance of an avalanche.")<jupyter_output>A weak_layer with value 5 yields a 56.51% chance of an avalanche.
A weak_layer with value 4 yields a 50.95% chance of an avalanche.
A weak_layer with value 7 yields a 67.05% chance of an avalanche.
A weak_layer with value 0 yields a 29.76% chance of an avalanche.<jupyter_text>Let's plot out model to understand this:<jupyter_code># plot the model
# Show a graph of the result
predict = lambda x: model.predict(pandas.DataFrame({"weak_layers": x}))

graphing.line_2D([("Model", predict)],
                 x_range=[-20,40],
                 label_x="weak_layers", 
                 label_y="estimated probability of an avalanche")<jupyter_output><empty_output><jupyter_text>The line plots the function of the __probability__ of an avalanche over the number of weak layers; Notice that the more weak layers, the more likely an avalanche will happen. This plot can look a bit confusing for two reasons.Firstly, the curve can make predictions from negative to positive infinity, but we only have data for 0 - 10 layers:<jupyter_code>print("Minimum number of weak layers:", min(train.weak_layers))
print("Maximum number of weak layers:", max(train.weak_layers))<jupyter_output>Minimum number of weak layers: 0
Maximum number of weak layers: 10<jupyter_text>This is because logistic regression models allow predictions outside the range of values they have seen, and sometimes do so quite well.The second reason the plot is confusing is that at 0 layers, there's still _some_ risk of an avalanche. Similarly, at 10 layers, there isn't a 100% risk of an avalanche. This is actually in line with the data:<jupyter_code>import numpy as np

# Get actual rates of avalanches at 0 years
avalanche_outcomes_for_0_layers = train[train.weak_layers == 0].avalanche
print("Average rate of avalanches for 0 weak layers of snow", np.average(avalanche_outcomes_for_0_layers))

# Get actual rates of avalanches at 10 years
avalanche_outcomes_for_10_layers = train[train.weak_layers == 10].avalanche
print("Average rate of avalanches for 10 weak layers of snow", np.average(avalanche_outcomes_for_10_layers))<jupyter_output>Average rate of avalanches for 0 weak layers of snow 0.3880597014925373
Average rate of avalanches for 10 weak layers of snow 0.7761194029850746<jupyter_text>Our model is actually doing a good job! It's just that avalanches aren't _only_ caused by weak layers of snow. If we want to do better, we probably need to think about including other information in the model.  Classification or decision thresholdsTo return a binary category (`True` = "avalanche", `False` = "no avalanche") we need to define a *Classification Threshold* value. Any probability above that threshold is returned as the positive category, whereas values below it will be returned as the negative category.Let's see what happens if set our threshold to `0.5` (meaning that our model will return `True` whenever it calculates a chance above 50% of an avalanche happening):<jupyter_code># threshold to get an absolute value
threshold = 0.5

# Add classification to the samples we used before
for sample, pred in list(zip(samples,estimated_probabilities)):
    print(f"A weak_layer with value {sample} yields a chance of {pred * 100:.2f}% of an avalanche. Classification = {pred > threshold}")<jupyter_output>A weak_layer with value 5 yields a chance of 56.51% of an avalanche. Classification = True
A weak_layer with value 4 yields a chance of 50.95% of an avalanche. Classification = True
A weak_layer with value 7 yields a chance of 67.05% of an avalanche. Classification = True
A weak_layer with value 0 yields a chance of 29.76% of an avalanche. Classification = False<jupyter_text>Note that a `0.5` threshold is just a starting point that needs to be tuned depending on the data we're trying to classify. Performance on test setNow let's use our `test` dataset to perform a quick evaluation on how the model did. For now, we'll just look at how often we correctly predicted if there would be an avalanche or not<jupyter_code># Classify the mdel predictions using the threshold
predictions = model.predict(test) > threshold

# Compare the predictions to the actual outcomes in the dataset
accuracy = np.average(predictions == test.avalanche)

# Print the evaluation
print(f"The model correctly predicted outcomes {accuracy * 100:.2f}% of time.")<jupyter_output>The model correctly predicted outcomes 65.05% of time.
{'fraction_non_alphanumeric': 0.0937179128162046, 'fraction_numerical': 0.028283394007280874, 'mean_word_length': 4.2314453125, 'pattern_counts': {'":': 1, '<': 36, '<?xml version=': 0, '>': 49, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '575326', 'license': 'CC0-1.0', 'n_tokens_mistral': 3223, 'n_tokens_neox': 2907, 'n_words': 1217, 'path': '7-3-exercise-build-regression-model.ipynb', 'repo_name': 'Sooner0931/MS_ML_NBs', 'size': '10713'}