Part 0A: Neural Networks & The Learning Mechanism

# Traditional programming: We write the rules
def approve_loan(income, credit_score, employment_years):
    if credit_score < 600:
        return False  # Rule 1: Bad credit = reject
    if income < 30000 and employment_years < 2:
        return False  # Rule 2: Low income + short employment = reject
    if income > 100000:
        return True   # Rule 3: High income = approve
    # ... hundreds more rules

Traditional: Human writes rules → Computer follows rules → Output
ML:          Human provides examples → Computer learns rules → Output

# Machine learning paradigm: We provide examples, computer learns rules
training_data = [
    # (income, credit_score, employment_years) → approved?
    (75000, 720, 5, True),   # This person was approved
    (25000, 580, 1, False),  # This person was rejected
    (120000, 690, 3, True),  # This person was approved
    (45000, 650, 8, True),   # This person was approved
    # ... thousands more examples
]
model = learn_from_examples(training_data)  # Computer figures out the rules
result = model.predict(55000, 670, 4)       # Apply learned rules to new case

Income contribution:      45000 × 0.00001  =  0.45
Credit contribution:      680 × 0.01       =  6.80
Employment contribution:  3 × 0.1          =  0.30
Debt contribution:        12000 × -0.00005 = -0.60
Bias:                                      = -3.00
─────────────────────────────────────────────────
Total (z):                                 =  3.95

import numpy as np

def single_neuron_step_by_step(inputs, weights, bias):
    """
    A single artificial neuron - the fundamental building block of neural networks.

    This function mimics how a loan officer might make a decision:
                1. Consider multiple factors (inputs)
                2. Weight each factor by its importance (weights)
                3. Combine into a single score (weighted sum + bias)
                4. Make a decision based on the score (activation)

    Parameters:
                -----------
    inputs : list or numpy array
        The raw data we're making a decision about.
        Example: [income, credit_score, years_employed, existing_debt]

    weights : list or numpy array
        How important each input is. Learned during training.
        Positive weight = input helps, Negative weight = input hurts
        Example: [0.00001, 0.01, 0.1, -0.00005]

    bias : float
        A baseline adjustment. Like a default tendency toward yes or no.
        Positive bias = lean toward yes, Negative bias = lean toward no
        Example: -3.0 (slightly conservative)

    Returns:
                --------
    tuple: (weighted_sum, pre_activation, output)
                                - weighted_sum: The sum of (input × weight) for each input
                                - pre_activation: weighted_sum + bias (the raw "score")
                                - output: The final decision (1 or 0)
    """

    # =========================================================================
    # STEP 1: Calculate each input's contribution
    # =========================================================================
    # Each input contributes to the decision based on its weight.
    # Think: "How much does this factor push toward yes or no?"

    contributions = []
    print("Step 1: Calculate each input's contribution")
    print("-" * 50)

    for i in range(len(inputs)):
        # Multiply input by its weight to get its contribution
        contribution = inputs[i] * weights[i]
        contributions.append(contribution)

        # Explain what's happening
        if weights[i] > 0:
            direction = "pushes toward YES"
        else:
            direction = "pushes toward NO"

        print(f"  Input {i+1}: {inputs[i]:>10.2f} × weight {weights[i]:>10.5f} = {contribution:>8.2f} ({direction})")

    # =========================================================================
    # STEP 2: Sum all contributions (weighted sum)
    # =========================================================================
    # Add up all the individual contributions to get a combined "evidence" score

    weighted_sum = sum(contributions)
    print(f"\nStep 2: Sum all contributions")
    print("-" * 50)
    print(f"  Weighted sum = {weighted_sum:.2f}")

    # =========================================================================
    # STEP 3: Add bias (pre-activation)
    # =========================================================================
    # The bias shifts our decision threshold.
    # Positive bias = easier to say yes, Negative bias = harder to say yes

    pre_activation = weighted_sum + bias
    print(f"\nStep 3: Add bias")
    print("-" * 50)
    print(f"  Pre-activation = weighted_sum + bias")
    print(f"  Pre-activation = {weighted_sum:.2f} + ({bias:.2f}) = {pre_activation:.2f}")

    # =========================================================================
    # STEP 4: Apply activation function (make decision)
    # =========================================================================
    # Convert the raw score into a final decision.
    # Here we use a simple "step function": positive = yes, negative/zero = no

    if pre_activation > 0:
        output = 1  # Approve
        decision = "APPROVE (score > 0)"
    else:
        output = 0  # Reject
        decision = "REJECT (score ≤ 0)"

    print(f"\nStep 4: Apply activation (decision rule)")
    print("-" * 50)
    print(f"  Rule: If pre-activation > 0, output 1 (yes), else output 0 (no)")
    print(f"  {pre_activation:.2f} > 0? → {decision}")

    return weighted_sum, pre_activation, output


# =============================================================================
# EXAMPLE: Loan approval decision
# =============================================================================

print("=" * 60)
print("LOAN APPROVAL NEURON EXAMPLE")
print("=" * 60)
print()

# Our applicant's data
inputs = [
    45000,   # Annual income (€)
    680,     # Credit score (300-850 scale)
    3,       # Years employed
    12000    # Existing debt (€)
]

# How important each factor is (these would be LEARNED in a real model)
weights = [
                0.00001,   # Income: small weight because numbers are large
                0.01,      # Credit score: moderate importance
                0.1,       # Years employed: some importance
                -0.00005   # Existing debt: NEGATIVE because debt is bad
]

# Baseline tendency (slightly conservative bank)
bias = -3.0

print("Applicant data:")
print(f"  Income: €{inputs[0]:,}")
print(f"  Credit Score: {inputs[1]}")
print(f"  Years Employed: {inputs[2]}")
print(f"  Existing Debt: €{inputs[3]:,}")
print()

# Run the neuron
weighted_sum, pre_activation, output = single_neuron_step_by_step(inputs, weights, bias)

print()
print("=" * 60)
print(f"FINAL DECISION: {'APPROVED' if output == 1 else 'REJECTED'}")
print("=" * 60)

============================================================
LOAN APPROVAL NEURON EXAMPLE
============================================================

Applicant data:
  Income: €45,000
  Credit Score: 680
  Years Employed: 3
  Existing Debt: €12,000

Step 1: Calculate each input's contribution
--------------------------------------------------
  Input 1:   45000.00 × weight    0.00001 =     0.45 (pushes toward YES)
  Input 2:     680.00 × weight    0.01000 =     6.80 (pushes toward YES)
  Input 3:       3.00 × weight    0.10000 =     0.30 (pushes toward YES)
  Input 4:   12000.00 × weight   -0.00005 =    -0.60 (pushes toward NO)

Step 2: Sum all contributions
--------------------------------------------------
  Weighted sum = 6.95

Step 3: Add bias
--------------------------------------------------
  Pre-activation = weighted_sum + bias
  Pre-activation = 6.95 + (-3.00) = 3.95

Step 4: Apply activation (decision rule)
--------------------------------------------------
  Rule: If pre-activation > 0, output 1 (yes), else output 0 (no)
        3.95 > 0? → APPROVE (score > 0)

============================================================
FINAL DECISION: APPROVED
============================================================

import numpy as np

def neuron(x, w, b, activation_fn):
    """
    A single neuron in compact form.

    Parameters:
                -----------
    x : numpy array, shape (n_features,)
        Input values
    w : numpy array, shape (n_features,)
        Weights (learned parameters)
    b : float
        Bias (learned parameter)
    activation_fn : function
        The activation function to apply

    Returns:
                --------
    float: The neuron's output
    """
    # np.dot computes the dot product (sum of element-wise products)
    # This is the "weighted sum" step in one line
    z = np.dot(w, x) + b

    # Apply activation function to get final output
    y = activation_fn(z)

    return y


# Simple step activation function
def step_activation(z):
    """Return 1 if z > 0, else 0"""
    return 1 if z > 0 else 0


# Example usage
x = np.array([45000, 680, 3, 12000])
w = np.array([0.00001, 0.01, 0.1, -0.00005])
b = -3.0

output = neuron(x, w, b, step_activation)
print(f"Output: {output}")  # Output: 1 (approved)

  ---

Neuron 1: z₁ = w₁₁x₁ + w₁₂x₂ + b₁
Neuron 2 (takes z₁ as input): z₂ = w₂₁z₁ + b₂

z₂ = w₂₁(w₁₁x₁ + w₁₂x₂ + b₁) + b₂
   = w₂₁w₁₁x₁ + w₂₁w₁₂x₂ + w₂₁b₁ + b₂
   = (w₂₁w₁₁)x₁ + (w₂₁w₁₂)x₂ + (w₂₁b₁ + b₂)
   = W₁x₁ + W₂x₂ + B

Neuron 1: z₁ = w₁₁x₁ + w₁₂x₂ + b₁
          a₁ = σ(z₁)           ← Non-linear activation!

Neuron 2: z₂ = w₂₁a₁ + b₂      ← Now takes σ(z₁), not z₁
          a₂ = σ(z₂)

import numpy as np

def sigmoid(z):
    """
    Sigmoid activation function.

    What it does:
                - Takes any real number as input
                - Outputs a value between 0 and 1
                - Smooth, S-shaped curve

    When to use:
                - Output layer for binary classification (probability of class 1)
                - Historically used in hidden layers, but rarely now (see ReLU)

    Why it works:
                - Large negative z → output ≈ 0
                - Large positive z → output ≈ 1
                - z = 0 → output = 0.5

    Problem:
                - For very large or very small z, the gradient is nearly 0
                - This causes "vanishing gradients" in deep networks (explained later)
    """
    return 1 / (1 + np.exp(-z))

# Examples
print("Sigmoid examples:")
print(f"  sigmoid(-10) = {sigmoid(-10):.6f}")  # ≈ 0 (very negative → near 0)
print(f"  sigmoid(-2)  = {sigmoid(-2):.6f}")   # = 0.119
print(f"  sigmoid(0)   = {sigmoid(0):.6f}")    # = 0.5 (midpoint)
print(f"  sigmoid(2)   = {sigmoid(2):.6f}")    # = 0.881
print(f"  sigmoid(10)  = {sigmoid(10):.6f}")   # ≈ 1 (very positive → near 1)

def relu(z):
    """
    ReLU (Rectified Linear Unit) activation function.

    What it does:
                - If z > 0, output z (pass through unchanged)
                - If z ≤ 0, output 0 (block negative values)

    When to use:
                - Default choice for hidden layers in modern networks
                - Works well in deep networks

    Why it's great:
                - Computationally efficient (just a comparison)
                - Doesn't saturate for positive values (gradient = 1)
                - Creates sparse representations (many zeros)

    Problem:
                - "Dead neurons": If a neuron's output is always negative,
      it will always output 0 and never update (gradient = 0)
    """
    return np.maximum(0, z)

# Examples
print("\nReLU examples:")
print(f"  relu(-10) = {relu(-10)}")   # 0 (negative → blocked)
print(f"  relu(-2)  = {relu(-2)}")    # 0
print(f"  relu(0)   = {relu(0)}")     # 0
print(f"  relu(2)   = {relu(2)}")     # 2 (positive → passed through)
print(f"  relu(10)  = {relu(10)}")    # 10

def softmax(z):
    """
    Softmax activation function.

    What it does:
                - Takes a vector of K raw scores (one per class)
                - Converts to K probabilities that sum to 1
                - Higher scores get higher probabilities

    When to use:
                - Output layer for multi-class classification
                - Example: Classifying images into 10 digit classes (0-9)

    Why it works:
                - Exponential makes all values positive
                - Division by sum ensures they add to 1
                - Preserves ranking (highest score → highest probability)

    Note: We subtract max(z) for numerical stability to prevent
    overflow when computing exp() of large numbers.
    """
    # Subtract max for numerical stability (doesn't change result)
    exp_z = np.exp(z - np.max(z))
    return exp_z / np.sum(exp_z)

# Example: Classifying a digit image
# Raw scores from a network for classes 0-9
scores = np.array([1.2, 0.5, 0.1, 3.8, 0.2, 0.1, 0.3, 0.9, 0.4, 0.2])

probabilities = softmax(scores)
print("\nSoftmax example (digit classification):")
print("Raw scores:", scores)
print("Probabilities:", np.round(probabilities, 3))
print(f"Sum of probabilities: {np.sum(probabilities):.6f}")  # Should be 1.0
print(f"Predicted class: {np.argmax(probabilities)} (highest probability)")

Single neuron decision boundary:

    Credit Score
         ↑
         |      ✓ ✓ ✓ ✓
         |    ✓ ✓ ✓ ✓ ✓
         |  ─────────────── ← Linear boundary
         |  ✗ ✗ ✗ ✗
         |✗ ✗ ✗ ✗ ✗
         └──────────────→ Income

Real fraud detection patterns:

    Transaction Amount
         ↑
         |    ✓   ✗ ✗   ✓
         |  ✓ ✓ ✗ ✗ ✗ ✓ ✓
         |    ✓ ✗ ✗ ✗ ✓    ← Can't separate with a line!
         |      ✓ ✓ ✓
         |        ✓
         └──────────────→ Transaction Time

import numpy as np

class SimpleNetwork:
    """
    A simple 2-layer neural network for binary classification.

    Architecture: Input → Hidden Layer (ReLU) → Output Layer (Sigmoid)

    This network can learn non-linear decision boundaries because:
                1. The hidden layer creates multiple non-linear features
                2. The output layer combines them for the final decision
    """

    def __init__(self, input_size, hidden_size, output_size):
        """
        Initialize the network with random weights.

        Parameters:
                                -----------
        input_size : int
            Number of input features (e.g., 2 for [income, credit_score])
        hidden_size : int
            Number of neurons in hidden layer (more = more capacity)
        output_size : int
            Number of output neurons (1 for binary classification)
        """
        # =====================================================================
        # Weight Initialization
        # =====================================================================
        # Weights connect neurons between layers.
        # Shape: (neurons_in_this_layer, neurons_in_previous_layer)
        #
        # We initialize with small random values. The scaling factor
        # (np.sqrt(2/n)) is called "He initialization" and helps training.
        # We'll explain WHY this matters in the training section.

        # Hidden layer weights: connects input to hidden
        # Shape: (hidden_size, input_size) = (4, 2) if hidden=4, input=2
        self.W1 = np.random.randn(hidden_size, input_size) * np.sqrt(2.0 / input_size)

        # Hidden layer biases: one per hidden neuron
        # Shape: (hidden_size, 1) = (4, 1)
        self.b1 = np.zeros((hidden_size, 1))

        # Output layer weights: connects hidden to output
        # Shape: (output_size, hidden_size) = (1, 4)
        self.W2 = np.random.randn(output_size, hidden_size) * np.sqrt(2.0 / hidden_size)

        # Output layer biases: one per output neuron
        # Shape: (output_size, 1) = (1, 1)
        self.b2 = np.zeros((output_size, 1))

        print(f"Network initialized:")
        print(f"  Input size: {input_size}")
        print(f"  Hidden size: {hidden_size}")
        print(f"  Output size: {output_size}")
        print(f"  W1 shape: {self.W1.shape} (hidden × input)")
        print(f"  b1 shape: {self.b1.shape}")
        print(f"  W2 shape: {self.W2.shape} (output × hidden)")
        print(f"  b2 shape: {self.b2.shape}")
        print(f"  Total parameters: {self.W1.size + self.b1.size + self.W2.size + self.b2.size}")

    def forward(self, X, verbose=False):
        """
        Forward pass: compute output given input.

        Parameters:
                                -----------
        X : numpy array, shape (input_size, batch_size)
            Input data. Each column is one example.
            Example: shape (2, 32) = 2 features, 32 examples

        verbose : bool
            If True, print intermediate steps

        Returns:
                                --------
        numpy array, shape (output_size, batch_size)
            Output predictions (probabilities for each example)
        """
        if verbose:
            print("\n" + "="*60)
            print("FORWARD PASS")
            print("="*60)
            print(f"\nInput X shape: {X.shape}")
            print(f"Input X:\n{X}")

        # =====================================================================
        # LAYER 1: Input → Hidden
        # =====================================================================

        # Step 1a: Compute weighted sum for hidden layer
        # Z1 = W1 @ X + b1
        # Matrix multiplication: (hidden, input) @ (input, batch) = (hidden, batch)
        # This computes ALL hidden neurons for ALL examples at once!
        Z1 = np.dot(self.W1, X) + self.b1

        if verbose:
            print(f"\n--- Layer 1 (Input → Hidden) ---")
            print(f"W1 @ X + b1 = Z1")
            print(f"({self.W1.shape}) @ ({X.shape}) + ({self.b1.shape}) = ({Z1.shape})")
            print(f"Z1 (pre-activation):\n{Z1}")

        # Step 1b: Apply ReLU activation
        # ReLU(z) = max(0, z) for each element
        A1 = np.maximum(0, Z1)

        if verbose:
            print(f"A1 = ReLU(Z1):\n{A1}")

        # =====================================================================
        # LAYER 2: Hidden → Output
        # =====================================================================

        # Step 2a: Compute weighted sum for output layer
        Z2 = np.dot(self.W2, A1) + self.b2

        if verbose:
            print(f"\n--- Layer 2 (Hidden → Output) ---")
            print(f"W2 @ A1 + b2 = Z2")
            print(f"({self.W2.shape}) @ ({A1.shape}) + ({self.b2.shape}) = ({Z2.shape})")
            print(f"Z2 (pre-activation):\n{Z2}")

        # Step 2b: Apply sigmoid activation
        # Sigmoid squashes to (0, 1) for probability output
        A2 = 1 / (1 + np.exp(-Z2))

        if verbose:
            print(f"A2 = sigmoid(Z2):\n{A2}")
            print(f"\nFinal output shape: {A2.shape}")
            print(f"Predictions (probabilities): {A2.flatten()}")

        # Store intermediate values (needed for backpropagation later)
        self.cache = {
            'X': X, 'Z1': Z1, 'A1': A1, 'Z2': Z2, 'A2': A2
        }

        return A2


# =============================================================================
# EXAMPLE: Forward pass demonstration
# =============================================================================

print("="*70)
print("FORWARD PASS DEMONSTRATION")
print("="*70)

# Create a small network
np.random.seed(42)  # For reproducible results
net = SimpleNetwork(input_size=2, hidden_size=4, output_size=1)

# Create sample input: 3 loan applications
# Each column is one application: [income (scaled), credit_score (scaled)]
X = np.array([
    [0.45, 0.80, 0.25],   # Incomes (scaled to ~0-1)
    [0.68, 0.90, 0.55]    # Credit scores (scaled to ~0-1)
])

print("\n" + "="*70)
print("Processing 3 loan applications...")
print("="*70)

# Run forward pass with detailed output
predictions = net.forward(X, verbose=True)

print("\n" + "="*70)
print("INTERPRETATION")
print("="*70)
print(f"Application 1 (income=0.45, credit=0.68): {predictions[0,0]:.1%} approval probability")
print(f"Application 2 (income=0.80, credit=0.90): {predictions[0,1]:.1%} approval probability")
print(f"Application 3 (income=0.25, credit=0.55): {predictions[0,2]:.1%} approval probability")

======================================================================
FORWARD PASS DEMONSTRATION
======================================================================
Network initialized:
  Input size: 2
  Hidden size: 4
  Output size: 1
  W1 shape: (4, 2) (hidden × input)
  b1 shape: (4, 1)
  W2 shape: (1, 4) (output × hidden)
  b2 shape: (1, 1)
  Total parameters: 17

======================================================================
Processing 3 loan applications...
======================================================================

============================================================
FORWARD PASS
============================================================

Input X shape: (2, 3)
Input X:
[[0.45 0.8  0.25]
 [0.68 0.9  0.55]]
    --- Layer 1 (Input → Hidden) ---
W1 @ X + b1 = Z1
((4, 2)) @ ((2, 3)) + ((4, 1)) = ((4, 3))
Z1 (pre-activation):
[[ 0.12950164  0.27293345  0.04813317]
 [ 1.32712014  1.8888777   0.99958856]
 [-0.26458215 -0.39804596 -0.18731367]
 [ 1.23250138  1.95406151  0.8168923 ]]
A1 = ReLU(Z1):
[[0.12950164 0.27293345 0.04813317]
 [1.32712014 1.8888777  0.99958856]
 [0.         0.         0.        ]
 [1.23250138 1.95406151 0.8168923 ]]
    --- Layer 2 (Hidden → Output) ---
W2 @ A1 + b2 = Z2
((1, 4)) @ ((4, 3)) + ((1, 1)) = ((1, 3))
Z2 (pre-activation):
[[ 0.06026819 -0.00945422  0.09849182]]
A2 = sigmoid(Z2):
[[0.51506249 0.49763646 0.52460307]]

Final output shape: (1, 3)
Predictions (probabilities): [0.51506249 0.49763646 0.52460307]

======================================================================
INTERPRETATION
======================================================================
Application 1 (income=0.45, credit=0.68): 51.5% approval probability
Application 2 (income=0.80, credit=0.90): 49.8% approval probability
Application 3 (income=0.25, credit=0.55): 52.5% approval probability

Instead of:
    for each example:
        for each neuron:
            compute weighted sum

We do:
    Z = W @ X + b  (one matrix operation, processes everything)

def mean_squared_error(y_true, y_pred):
    """
    Mean Squared Error loss for regression problems.

    How it works:
                1. For each prediction, compute error: (actual - predicted)
                2. Square each error (makes all errors positive, penalizes big errors more)
                3. Average across all examples

    Parameters:
                -----------
    y_true : numpy array
        Actual values (ground truth)
    y_pred : numpy array
        Model's predictions

    Returns:
                --------
    float : The average squared error

    Example:
                --------
    Predicting house prices (in $100,000s):
    Actual:    [3.0, 4.5, 2.0]  ($300K, $450K, $200K)
    Predicted: [2.8, 4.2, 2.5]  ($280K, $420K, $250K)

    Errors:    [0.2, 0.3, -0.5]
    Squared:   [0.04, 0.09, 0.25]
    MSE:       0.1267

    Interpretation: On average, our squared error is 0.127, meaning
    typical error is roughly sqrt(0.127) ≈ 0.36 ($36,000)
    """
    # Step 1: Compute errors
    errors = y_true - y_pred

    # Step 2: Square the errors
    # This makes negative errors positive AND penalizes large errors more
    squared_errors = errors ** 2

    # Step 3: Average across all examples
    mse = np.mean(squared_errors)

    return mse


# Example: House price prediction
y_true = np.array([3.0, 4.5, 2.0])  # Actual prices
y_pred = np.array([2.8, 4.2, 2.5])  # Predicted prices

mse = mean_squared_error(y_true, y_pred)
print(f"MSE: {mse:.4f}")
print(f"Typical error magnitude: ${np.sqrt(mse) * 100000:,.0f}")

def binary_cross_entropy(y_true, y_pred, epsilon=1e-15):
    """
    Binary Cross-Entropy loss for binary classification.

    How it works:
                - If true label is 1, we want predicted probability to be HIGH
      Loss = -log(predicted). If predicted=0.9, loss=-log(0.9)=0.1 (low, good!)
                              If predicted=0.1, loss=-log(0.1)=2.3 (high, bad!)
                    - If true label is 0, we want predicted probability to be LOW
      Loss = -log(1-predicted). If predicted=0.1, loss=-log(0.9)=0.1 (low, good!)
                                If predicted=0.9, loss=-log(0.1)=2.3 (high, bad!)

    Key insight: The -log function heavily penalizes confident WRONG predictions.
                - Predicting 0.99 when true label is 0 → loss ≈ 4.6 (very bad!)
                - Predicting 0.51 when true label is 0 → loss ≈ 0.7 (not great, but not terrible)

    This encourages the model to only be confident when it's right.

    Parameters:
                -----------
    y_true : numpy array
        Actual labels (0 or 1)
    y_pred : numpy array
        Predicted probabilities (between 0 and 1)
    epsilon : float
        Small value to prevent log(0) which is undefined

    Returns:
                --------
    float : The average binary cross-entropy loss
    """
    # Clip predictions to prevent log(0) or log(1)
    # log(0) = -infinity, which would break our computation
    y_pred = np.clip(y_pred, epsilon, 1 - epsilon)

    # Compute loss for each example
    # When y_true=1: loss = -log(y_pred)
    # When y_true=0: loss = -log(1 - y_pred)
    # This formula handles both cases elegantly:
    individual_losses = -(
        y_true * np.log(y_pred) +           # Active when y_true=1
        (1 - y_true) * np.log(1 - y_pred)   # Active when y_true=0
    )

    # Average across all examples
    bce = np.mean(individual_losses)

    return bce


# Example: Loan approval classification
print("Binary Cross-Entropy Examples:")
print("-" * 50)

# Case 1: Correct and confident prediction
y_true = np.array([1])       # Should be approved
y_pred = np.array([0.95])    # Model is confident: 95% approval
loss = binary_cross_entropy(y_true, y_pred)
print(f"True=1, Pred=0.95 (confident, correct): Loss = {loss:.4f}")

# Case 2: Correct but uncertain prediction
y_pred = np.array([0.55])    # Model is uncertain: 55% approval
loss = binary_cross_entropy(y_true, y_pred)
print(f"True=1, Pred=0.55 (uncertain, correct): Loss = {loss:.4f}")

# Case 3: Wrong and confident prediction (VERY BAD)
y_pred = np.array([0.05])    # Model is confident it should be rejected!
loss = binary_cross_entropy(y_true, y_pred)
print(f"True=1, Pred=0.05 (confident, WRONG): Loss = {loss:.4f}")

print()
print("Notice: Confident wrong predictions have MUCH higher loss!")
print("This teaches the model to be uncertain when it doesn't know.")

Input → [W1] → Z1 → [ReLU] → A1 → [W2] → Z2 → [Sigmoid] → A2 → [Loss] → L

import numpy as np

def backprop_step_by_step():
    """
    Demonstrate backpropagation with concrete numbers.

    Network: 2 inputs → 2 hidden neurons → 1 output
    """
    np.random.seed(42)

    print("="*70)
    print("BACKPROPAGATION STEP BY STEP")
    print("="*70)

    # =========================================================================
    # SETUP: Define a simple network with known weights
    # =========================================================================

    # Network architecture: 2 → 2 → 1

    # Hidden layer weights (2 neurons, each receiving 2 inputs)
    W1 = np.array([
        [0.1, 0.2],   # Neuron 1 weights
        [0.3, 0.4]    # Neuron 2 weights
    ])
    b1 = np.array([[0.1], [0.2]])  # Biases for hidden neurons

    # Output layer weights (1 neuron receiving 2 inputs from hidden)
    W2 = np.array([[0.5, 0.6]])  # Single output neuron
    b2 = np.array([[0.3]])       # Bias for output

    print("\n--- NETWORK PARAMETERS ---")
    print(f"W1 (hidden layer weights):\n{W1}")
    print(f"b1 (hidden layer biases): {b1.flatten()}")
    print(f"W2 (output layer weights): {W2}")
    print(f"b2 (output layer bias): {b2.flatten()}")

    # =========================================================================
    # INPUT DATA
    # =========================================================================

    # Single training example
    X = np.array([[0.5], [0.8]])  # Input features
    Y = np.array([[1.0]])          # True label (should output 1)

    print(f"\n--- INPUT ---")
    print(f"X (input features): {X.flatten()}")
    print(f"Y (true label): {Y.flatten()}")

    # =========================================================================
    # FORWARD PASS (compute predictions)
    # =========================================================================

    print("\n" + "="*70)
    print("FORWARD PASS")
    print("="*70)

    # --- Hidden Layer ---
    # Z1 = W1 @ X + b1
    Z1 = np.dot(W1, X) + b1
    print(f"\n[Hidden Layer] Z1 = W1 @ X + b1")
    print(f"  W1 @ X = {W1} @ {X.flatten()} = {np.dot(W1, X).flatten()}")
    print(f"  Z1 = {np.dot(W1, X).flatten()} + {b1.flatten()} = {Z1.flatten()}")

    # Apply ReLU: A1 = max(0, Z1)
    A1 = np.maximum(0, Z1)
    print(f"\n[Hidden Layer] A1 = ReLU(Z1) = max(0, Z1)")
    print(f"  A1 = max(0, {Z1.flatten()}) = {A1.flatten()}")

    # --- Output Layer ---
    # Z2 = W2 @ A1 + b2
    Z2 = np.dot(W2, A1) + b2
    print(f"\n[Output Layer] Z2 = W2 @ A1 + b2")
    print(f"  Z2 = {W2} @ {A1.flatten()} + {b2.flatten()} = {Z2.flatten()}")

    # Apply Sigmoid: A2 = 1 / (1 + exp(-Z2))
    A2 = 1 / (1 + np.exp(-Z2))
    print(f"\n[Output Layer] A2 = sigmoid(Z2)")
    print(f"  A2 = sigmoid({Z2.flatten()}) = {A2.flatten()}")

    # --- Compute Loss ---
    # Binary Cross-Entropy: L = -[Y*log(A2) + (1-Y)*log(1-A2)]
    epsilon = 1e-15
    L = -np.mean(Y * np.log(A2 + epsilon) + (1 - Y) * np.log(1 - A2 + epsilon))
    print(f"\n[Loss] Binary Cross-Entropy")
    print(f"  L = {L:.6f}")
    print(f"  Interpretation: Predicted {A2[0,0]:.4f}, true label is {Y[0,0]}")
    print(f"  Error: We predicted {A2[0,0]:.1%} probability, should be {Y[0,0]:.0%}")

    # =========================================================================
    # BACKWARD PASS (compute gradients)
    # =========================================================================

    print("\n" + "="*70)
    print("BACKWARD PASS")
    print("="*70)
    print("\nWe now compute: How does changing each weight affect the loss?")

    m = 1  # Number of examples (batch size)

    # --- Step 1: Gradient of Loss w.r.t. A2 ---
    # For BCE: dL/dA2 = -Y/A2 + (1-Y)/(1-A2)
    dL_dA2 = -(Y / (A2 + epsilon)) + (1 - Y) / (1 - A2 + epsilon)
    print(f"\n[Step 1] Gradient of Loss w.r.t. Output (dL/dA2)")
    print(f"  Formula: dL/dA2 = -Y/A2 + (1-Y)/(1-A2)")
    print(f"  dL/dA2 = {dL_dA2.flatten()[0]:.6f}")
    print(f"  Meaning: If we increase A2 slightly, loss changes by this amount")

    # --- Step 2: Gradient through Sigmoid ---
    # Sigmoid derivative: dA2/dZ2 = A2 * (1 - A2)
    dA2_dZ2 = A2 * (1 - A2)
    dL_dZ2 = dL_dA2 * dA2_dZ2  # Chain rule!
    print(f"\n[Step 2] Chain through Sigmoid (dL/dZ2 = dL/dA2 × dA2/dZ2)")
    print(f"  Sigmoid derivative: dA2/dZ2 = A2 × (1-A2) = {A2[0,0]:.4f} × {1-A2[0,0]:.4f} = {dA2_dZ2[0,0]:.6f}")
    print(f"  Chain rule: dL/dZ2 = {dL_dA2[0,0]:.6f} × {dA2_dZ2[0,0]:.6f} = {dL_dZ2[0,0]:.6f}")

    # --- Step 3: Gradients for Output Layer Weights ---
    # dL/dW2 = dL/dZ2 × dZ2/dW2 = dL/dZ2 × A1
    dL_dW2 = np.dot(dL_dZ2, A1.T) / m
    dL_db2 = np.sum(dL_dZ2, axis=1, keepdims=True) / m
    print(f"\n[Step 3] Gradients for Output Layer (W2, b2)")
    print(f"  dL/dW2 = dL/dZ2 × A1ᵀ = {dL_dZ2[0,0]:.6f} × {A1.flatten()} = {dL_dW2.flatten()}")
    print(f"  dL/db2 = dL/dZ2 = {dL_db2.flatten()}")
    print(f"  Meaning: These tell us how to adjust W2 and b2 to reduce loss")

    # --- Step 4: Backpropagate to Hidden Layer ---
    # dL/dA1 = W2ᵀ × dL/dZ2
    dL_dA1 = np.dot(W2.T, dL_dZ2)
    print(f"\n[Step 4] Backpropagate error to Hidden Layer (dL/dA1)")
    print(f"  dL/dA1 = W2ᵀ × dL/dZ2 = {W2.T.flatten()} × {dL_dZ2[0,0]:.6f} = {dL_dA1.flatten()}")
    print(f"  Meaning: This is how much each hidden neuron contributed to the error")

    # --- Step 5: Gradient through ReLU ---
    # ReLU derivative: 1 if Z1 > 0, else 0
    dA1_dZ1 = (Z1 > 0).astype(float)  # 1 where Z1 > 0, else 0
    dL_dZ1 = dL_dA1 * dA1_dZ1  # Chain rule
    print(f"\n[Step 5] Chain through ReLU (dL/dZ1 = dL/dA1 × dA1/dZ1)")
    print(f"  ReLU derivative: dA1/dZ1 = {dA1_dZ1.flatten()} (1 if Z1>0, else 0)")
    print(f"  dL/dZ1 = {dL_dA1.flatten()} × {dA1_dZ1.flatten()} = {dL_dZ1.flatten()}")

    # --- Step 6: Gradients for Hidden Layer Weights ---
    dL_dW1 = np.dot(dL_dZ1, X.T) / m
    dL_db1 = np.sum(dL_dZ1, axis=1, keepdims=True) / m
    print(f"\n[Step 6] Gradients for Hidden Layer (W1, b1)")
    print(f"  dL/dW1 = dL/dZ1 × Xᵀ:")
    print(f"    {dL_dW1}")
    print(f"  dL/db1 = {dL_db1.flatten()}")

    # =========================================================================
    # WEIGHT UPDATE
    # =========================================================================

    print("\n" + "="*70)
    print("WEIGHT UPDATE (Gradient Descent)")
    print("="*70)

    learning_rate = 0.1
    print(f"\nLearning rate: {learning_rate}")
    print(f"Update rule: W_new = W_old - learning_rate × gradient")

    W2_new = W2 - learning_rate * dL_dW2
    b2_new = b2 - learning_rate * dL_db2
    W1_new = W1 - learning_rate * dL_dW1
    b1_new = b1 - learning_rate * dL_db1

    print(f"\nW2: {W2.flatten()} → {W2_new.flatten()}")
    print(f"b2: {b2.flatten()} → {b2_new.flatten()}")
    print(f"W1:\n{W1} →\n{W1_new}")
    print(f"b1: {b1.flatten()} → {b1_new.flatten()}")

    # =========================================================================
    # VERIFY IMPROVEMENT
    # =========================================================================

    print("\n" + "="*70)
    print("VERIFY: Did the update help?")
    print("="*70)

    # Forward pass with new weights
    Z1_new = np.dot(W1_new, X) + b1_new
    A1_new = np.maximum(0, Z1_new)
    Z2_new = np.dot(W2_new, A1_new) + b2_new
    A2_new = 1 / (1 + np.exp(-Z2_new))
    L_new = -np.mean(Y * np.log(A2_new + epsilon) + (1 - Y) * np.log(1 - A2_new + epsilon))

    print(f"\nBefore update:")
    print(f"  Prediction: {A2[0,0]:.4f} (target: {Y[0,0]:.1f})")
    print(f"  Loss: {L:.6f}")

    print(f"\nAfter update:")
    print(f"  Prediction: {A2_new[0,0]:.4f} (target: {Y[0,0]:.1f})")
    print(f"  Loss: {L_new:.6f}")

    print(f"\nImprovement: Loss decreased by {L - L_new:.6f} ({(L-L_new)/L*100:.1f}%)")
    print("The prediction moved closer to the target! 🎉")

    return {
        'W1': W1_new, 'b1': b1_new, 'W2': W2_new, 'b2': b2_new,
        'final_loss': L_new
    }


# Run the demonstration
result = backprop_step_by_step()

def learning_rate_demo():
    """
    Demonstrate the effect of learning rate.

    We'll minimize a simple function: f(x) = x²
    The minimum is at x = 0.
    """

    def f(x):
        return x ** 2

    def gradient(x):
        return 2 * x  # Derivative of x²

    print("="*60)
    print("LEARNING RATE COMPARISON")
    print("="*60)
    print("Minimizing f(x) = x², starting at x = 10")
    print("Optimal solution: x = 0")
    print()

    learning_rates = [0.01, 0.1, 0.5, 0.9, 5.0]

    for lr in learning_rates:
        x = 10.0  # Starting point
        history = [x]

        print(f"Learning rate = {lr}:")

        for step in range(20):
            grad = gradient(x)
            x = x - lr * grad  # Gradient descent update
            history.append(x)

            # Check for divergence
            if abs(x) > 1e10:
                print(f"  Step {step+1}: DIVERGED! (x = {x:.2e})")
                break
        else:
            print(f"  After 20 steps: x = {x:.6f}, f(x) = {f(x):.6f}")

        # Categorize behavior
        if abs(x) > 1e10:
            print(f"  Result: DIVERGED (learning rate too high)")
        elif abs(x) < 0.01:
            print(f"  Result: CONVERGED (good learning rate)")
        else:
            print(f"  Result: STILL CONVERGING (learning rate may be too low)")
        print()


learning_rate_demo()

class AdamOptimizer:
    """
    Adam optimizer - the practical default for most neural networks.

    Key ideas:
                1. Momentum: Keep a running average of past gradients.
       If gradients consistently point the same direction, speed up.
       This helps escape flat regions and smooth out noise.
                    2. Adaptive learning rates: Keep track of how much each parameter's
       gradient varies. Parameters with wildly varying gradients get
       smaller learning rates. Parameters with consistent gradients
       get larger learning rates.

    Default hyperparameters (rarely need changing):
                - learning_rate = 0.001
                - beta1 = 0.9 (momentum decay)
                - beta2 = 0.999 (squared gradient decay)
                - epsilon = 1e-8 (prevents division by zero)
    """

    def __init__(self, learning_rate=0.001, beta1=0.9, beta2=0.999, epsilon=1e-8):
        """
        Parameters:
                                -----------
        learning_rate : float
            Base step size. 0.001 is good for most problems.

        beta1 : float
            Decay rate for first moment (gradient average).
                                                0.9 means "remember 90% of past gradients".

        beta2 : float
            Decay rate for second moment (squared gradient average).
                                                0.999 means "remember 99.9% of past squared gradients".

        epsilon : float
            Small constant to prevent division by zero.
        """
        self.lr = learning_rate
        self.beta1 = beta1
        self.beta2 = beta2
        self.epsilon = epsilon

        # These store running averages for each parameter
        self.m = {}  # First moment (gradient average)
        self.v = {}  # Second moment (squared gradient average)
        self.t = 0   # Timestep counter

    def update(self, params, grads):
        """
        Update parameters using Adam algorithm.

        Parameters:
                                -----------
        params : dict
            Dictionary of parameters, e.g., {'W1': array, 'b1': array, ...}
        grads : dict
            Dictionary of gradients, e.g., {'dW1': array, 'db1': array, ...}
        """
        self.t += 1  # Increment timestep

        for key in params:
            # Get gradient for this parameter
            # Convert 'W1' to 'dW1', 'b1' to 'db1', etc.
            grad_key = 'd' + key
            if grad_key not in grads:
                continue

            grad = grads[grad_key]

            # Initialize moments on first update
            if key not in self.m:
                self.m[key] = np.zeros_like(params[key])
                self.v[key] = np.zeros_like(params[key])

            # ================================================================
            # Step 1: Update first moment (momentum)
            # ================================================================
            # m = β₁ × m + (1 - β₁) × gradient
            # This is a weighted average of past gradients
            self.m[key] = self.beta1 * self.m[key] + (1 - self.beta1) * grad

            # ================================================================
            # Step 2: Update second moment (adaptive learning rate)
            # ================================================================
            # v = β₂ × v + (1 - β₂) × gradient²
            # This tracks how much the gradient varies
            self.v[key] = self.beta2 * self.v[key] + (1 - self.beta2) * (grad ** 2)

            # ================================================================
            # Step 3: Bias correction
            # ================================================================
            # The moments are biased toward zero initially (they start at 0)
            # This correction compensates for that bias
            m_corrected = self.m[key] / (1 - self.beta1 ** self.t)**
**            v_corrected = self.v[key] / (1 - self.beta2 ** self.t)

            # ================================================================
            # Step 4: Update parameter
            # ================================================================
            # The key insight: we divide by sqrt(v), which gives each parameter
            # its own effective learning rate based on gradient history
            params[key] -= self.lr * m_corrected / (np.sqrt(v_corrected) + self.epsilon)

import numpy as np
def adam_demo():
    """
    Compare Adam vs vanilla gradient descent on f(x) = x²
    """
    def f(x):
        return x ** 2

    def gradient(x):
        return 2 * x

    print("="*60)
    print("ADAM vs VANILLA GRADIENT DESCENT")
    print("="*60)
    print("Minimizing f(x) = x², starting at x = 10\n")

    # Vanilla GD with a conservative learning rate
    x_gd = 10.0
    lr = 0.1
    gd_history = [x_gd]

    for _ in range(20):
        x_gd = x_gd - lr * gradient(x_gd)
        gd_history.append(x_gd)

    # Adam optimizer
    x_adam = 10.0
    adam = AdamOptimizer(learning_rate=0.5)  # Can use larger lr with Adam
    adam_history = [x_adam]

    for _ in range(20):
        # Adam expects dict format, we'll simulate it
        params = {'x': np.array([x_adam])}
        grads = {'dx': np.array([gradient(x_adam)])}
        adam.update(params, grads)
        x_adam = params['x'][0]
        adam_history.append(x_adam)

    print(f"Vanilla GD (lr={lr}):")
    print(f"  Step 5:  x = {gd_history[5]:.6f}")
    print(f"  Step 10: x = {gd_history[10]:.6f}")
    print(f"  Step 20: x = {gd_history[20]:.6f}")
    print()
    print(f"Adam (lr=0.5):")
    print(f"  Step 5:  x = {adam_history[5]:.6f}")
    print(f"  Step 10: x = {adam_history[10]:.6f}")
    print(f"  Step 20: x = {adam_history[20]:.6f}")
    print()
    print("Adam adapts its effective learning rate based on gradient history.")

adam_demo()

def demonstrate_vanishing_gradient():
    """
    Show why deep networks with sigmoid activations are hard to train.

    Sigmoid derivative: σ'(z) = σ(z)(1 - σ(z))
    Maximum value: 0.25 (when z = 0)

    In backpropagation, we multiply by this derivative at each layer.
    After many layers: 0.25 × 0.25 × 0.25 × ... → ~0
    """

    print("Gradient magnitude through layers with Sigmoid:")
    print("-" * 50)

    # Maximum sigmoid derivative
    sigmoid_deriv_max = 0.25

    gradient = 1.0  # Start with gradient = 1 from the loss

    for layer in range(1, 11):
        # Each layer multiplies gradient by sigmoid derivative
        gradient = gradient * sigmoid_deriv_max
        print(f"  Layer {layer:2d}: gradient magnitude = {gradient:.2e}")

    print()
    print(f"After 10 layers: gradient is {gradient:.2e}")
    print("This is essentially ZERO - early layers learn nothing!")
    print()
    print("This is why deep networks historically couldn't be trained,")
    print("until ReLU and other techniques came along.")


demonstrate_vanishing_gradient()