# ---------------------------------------------------------------------------
# penguin_species.py
#
# AI project on Statistical Data: predicting a penguin's species from four
# measurements of its body.
#
# CBSE Class X, Artificial Intelligence (417), Part C project work.
# SDG 14 (Life below Water) and SDG 15 (Life on Land).
#
# The program follows the six stages of the AI Project Cycle, in order:
#   1. Problem scoping      2. Data acquisition     3. Data exploration
#   4. Modelling            5. Evaluation           6. Deployment
#
# It builds two models and compares them:
#   - a RULE BASED model, whose rules we write ourselves after looking at the
#     charts, and
#   - a LEARNING BASED model, which works out the "typical" penguin of each
#     species from the training data on its own.
#
# Libraries: pandas (the data), numpy (the arithmetic), matplotlib (charts).
# ---------------------------------------------------------------------------

# %% Stage 1 - Problem scoping
print("=" * 60)
print("PENGUIN SPECIES PREDICTION")
print("=" * 60)
print("Who   : field scientists counting penguins in Antarctica")
print("What  : telling three similar species apart takes an expert")
print("Where : on the islands, during the breeding season")
print("Why   : species counts show how each one is coping with change")
print()

# %% Stage 2 - Data acquisition
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# Show every column when a table is printed, instead of "..." in the middle.
pd.set_option("display.width", 120)
pd.set_option("display.max_columns", 20)

data = pd.read_csv("penguins.csv")

print("--- Stage 2: Data acquisition ---")
print("Rows read    :", len(data))
print("Columns read :", len(data.columns))
print()
print("The first 10 rows:")
print(data.head(10))
print()
data.info()
print()

# %% Stage 3 - Data exploration
print("--- Stage 3: Data exploration ---")

# The four measurements the model will use. These are the FEATURES; the
# species is the LABEL - the answer the model has to predict.
features = ["bill_length_mm", "bill_depth_mm", "flipper_length_mm", "body_mass_g"]
label = "species"

# A penguin with a measurement missing cannot be compared with anything, so
# that row is removed. The count is printed so nothing disappears silently.
missing = data[features].isnull().any(axis=1).sum()
data = data.dropna(subset=features)
print("Rows with a measurement missing :", missing, "(removed)")
print("Rows left                       :", len(data))
print()

print("Penguins of each species:")
counts = data[label].value_counts()
print(counts)
print()

# Mean, median and mode of one feature, worked out with numpy and pandas.
flipper = data["flipper_length_mm"]
print("Flipper length - mean   :", round(np.mean(flipper), 1), "mm")
print("Flipper length - median :", np.median(flipper), "mm")
print("Flipper length - mode   :", flipper.mode()[0], "mm")
print()

print("Average measurements of each species:")
print(data.groupby(label)[features].mean().round(1))
print()

colours = {"Adelie": "#f28e2b", "Chinstrap": "#9467bd", "Gentoo": "#17a2b8"}
species_list = ["Adelie", "Chinstrap", "Gentoo"]

# Chart 1: how many of each species - is the data balanced?
plt.figure(figsize=(6, 4))
plt.bar(species_list, counts.reindex(species_list, fill_value=0),
        color=[colours["Adelie"], colours["Chinstrap"], colours["Gentoo"]])
plt.title("Penguins of each species in the dataset")
plt.xlabel("Species")
plt.ylabel("Number of penguins")
plt.tight_layout()
plt.savefig("chart1_species_count.png")
plt.close()

# Chart 2: flipper length against body mass, one colour per species.
plt.figure(figsize=(7, 5))
for s in species_list:
    rows = data[data[label] == s]
    plt.scatter(rows["flipper_length_mm"], rows["body_mass_g"],
                color=colours[s], label=s, alpha=0.75)
plt.title("Flipper length and body mass")
plt.xlabel("Flipper length (mm)")
plt.ylabel("Body mass (g)")
plt.legend()
plt.grid(True, linestyle="--", alpha=0.4)
plt.tight_layout()
plt.savefig("chart2_flipper_mass.png")
plt.close()

# Chart 3: bill length against bill depth.
plt.figure(figsize=(7, 5))
for s in species_list:
    rows = data[data[label] == s]
    plt.scatter(rows["bill_length_mm"], rows["bill_depth_mm"],
                color=colours[s], label=s, alpha=0.75)
plt.title("Bill length and bill depth")
plt.xlabel("Bill length (mm)")
plt.ylabel("Bill depth (mm)")
plt.legend()
plt.grid(True, linestyle="--", alpha=0.4)
plt.tight_layout()
plt.savefig("chart3_bill.png")
plt.close()

# %% Stage 4 - Modelling
print("--- Stage 4: Modelling ---")

# Train-test split. The rows are shuffled (random_state fixes the shuffle so
# every run is the same), then 80% are used to TRAIN and 20% are kept aside
# to TEST. The model never sees the testing rows while it is being built.
data = data.sample(frac=1, random_state=7)
cut = int(len(data) * 0.8)
train = data[:cut]
test = data[cut:]
print("Training rows :", len(train))
print("Testing rows  :", len(test))
print()

# Model A - RULE BASED. These two rules were written by hand after looking
# at charts 2 and 3: Gentoo have the longest flippers, and of the other two,
# Chinstrap have the longer bill.
rule_predictions = []
for i in range(len(test)):
    row = test.iloc[i]
    if row["flipper_length_mm"] >= 206:
        rule_predictions.append("Gentoo")
    elif row["bill_length_mm"] >= 44:
        rule_predictions.append("Chinstrap")
    else:
        rule_predictions.append("Adelie")

# Model B - LEARNING BASED. The model learns the average penguin of each
# species from the training rows, then labels a new penguin with whichever
# average it is closest to.
#
# Body mass is in thousands and bill depth is under 25, so measured as they
# are, mass would drown out everything else. Each feature is first rescaled
# to run from 0 to 1, using the smallest and largest value in TRAINING data.
low = train[features].min()
high = train[features].max()
train_scaled = (train[features] - low) / (high - low)
test_scaled = (test[features] - low) / (high - low)

averages = train_scaled.groupby(train[label]).mean()
print("What the model learned - the average penguin of each species")
print("(every feature rescaled to run from 0 to 1):")
print(averages.round(2))
print()

learned_predictions = []
for i in range(len(test)):
    row = test_scaled.iloc[i]
    # Distance from this penguin to each species' average: the gaps are
    # squared and added up. The smallest total is the closest species.
    distance = ((averages - row) ** 2).sum(axis=1)
    learned_predictions.append(distance.idxmin())

# %% Stage 5 - Evaluation
print("--- Stage 5: Evaluation ---")
actual = list(test[label])

rule_correct = 0
learned_correct = 0
for i in range(len(actual)):
    if rule_predictions[i] == actual[i]:
        rule_correct = rule_correct + 1
    if learned_predictions[i] == actual[i]:
        learned_correct = learned_correct + 1

rule_accuracy = rule_correct / len(actual) * 100
learned_accuracy = learned_correct / len(actual) * 100
print("Rule based model     :", rule_correct, "of", len(actual), "correct")
print("Learning based model :", learned_correct, "of", len(actual), "correct")
print("Accuracy of the rule based model     :", round(rule_accuracy, 1), "%")
print("Accuracy of the learning based model :", round(learned_accuracy, 1), "%")
print()

# The confusion matrix: rows are what the penguin really was, columns are
# what the model said. Everything on the diagonal is a correct prediction.
matrix = pd.crosstab(pd.Series(actual, name="Actual"),
                     pd.Series(learned_predictions, name="Predicted"))
matrix = matrix.reindex(index=species_list, columns=species_list, fill_value=0)
print("Confusion matrix of the learning based model:")
print(matrix)
print()

# Precision and recall for each species, read straight off the matrix.
print("Species      Precision   Recall")
for s in species_list:
    tp = matrix.loc[s, s]
    predicted_as_s = matrix[s].sum()
    really_s = matrix.loc[s].sum()
    if predicted_as_s > 0:
        precision = tp / predicted_as_s * 100
    else:
        precision = 0
    if really_s > 0:
        recall = tp / really_s * 100
    else:
        recall = 0
    print(s.ljust(12), str(round(precision, 1)).rjust(8), "%", str(round(recall, 1)).rjust(7), "%")
print()

# Chart 4: the confusion matrix, drawn as a grid of coloured boxes.
plt.figure(figsize=(5.5, 4.5))
plt.imshow(matrix.values, cmap="Blues")
for r in range(len(species_list)):
    for c in range(len(species_list)):
        value = matrix.values[r][c]
        if value > matrix.values.max() / 2:
            ink = "white"
        else:
            ink = "black"
        plt.text(c, r, value, ha="center", va="center", fontsize=14, color=ink)
plt.xticks(range(3), species_list)
plt.yticks(range(3), species_list)
plt.xlabel("Predicted species")
plt.ylabel("Actual species")
plt.title("Confusion matrix (testing data)")
plt.tight_layout()
plt.savefig("chart4_confusion_matrix.png")
plt.close()

# Chart 5: the two models side by side.
plt.figure(figsize=(6, 4))
plt.bar(["Rule based", "Learning based"], [rule_accuracy, learned_accuracy],
        color=["#999999", "#17a2b8"])
plt.ylim(0, 100)
plt.title("Accuracy on the same testing data")
plt.ylabel("Accuracy (%)")
plt.text(0, rule_accuracy + 2, str(round(rule_accuracy, 1)) + "%", ha="center")
plt.text(1, learned_accuracy + 2, str(round(learned_accuracy, 1)) + "%", ha="center")
plt.tight_layout()
plt.savefig("chart5_models_compared.png")
plt.close()

# %% Stage 6 - Deployment
print("--- Stage 6: Deployment ---")
print("Three penguins measured on a new visit:")

new_penguins = pd.DataFrame({
    "bill_length_mm":    [39.0, 50.5, 47.0],
    "bill_depth_mm":     [18.5, 19.0, 14.8],
    "flipper_length_mm": [188, 198, 215],
    "body_mass_g":       [3600, 3800, 5200],
})
new_scaled = (new_penguins - low) / (high - low)

for i in range(len(new_penguins)):
    distance = ((averages - new_scaled.iloc[i]) ** 2).sum(axis=1)
    m = new_penguins.iloc[i]
    print("Bill", m["bill_length_mm"], "x", m["bill_depth_mm"], "mm, flipper",
          int(m["flipper_length_mm"]), "mm,", int(m["body_mass_g"]), "g  ->",
          distance.idxmin())