Cross-Entropy

1
#include <bits/stdc++.h>
2
using namespace std;
3
4
// Normalize a vector to sum to 1 and clamp negatives to 0
5
vector<double> normalize(const vector<double>& v) {
6
    vector<double> out(v.size());
7
    double sum = 0.0;
8
    for (double x : v) sum += max(0.0, x);
9
    if (sum == 0.0) throw runtime_error("All entries are non-positive; cannot normalize.");
10
    for (size_t i = 0; i < v.size(); ++i) out[i] = max(0.0, v[i]) / sum;
11
    return out;
12
}
13
14
// Safe log with epsilon to avoid log(0).
15
double safe_log(double x, double eps = 1e-15) {
16
    if (x <= eps) x = eps;
17
    return log(x); // natural log -> nats
18
}
19
20
struct InfoMeasures {
21
    double entropy_P;     // H(P)
22
    double cross_entropy; // H(P,Q)
23
    double kl_PQ;         // D_KL(P||Q)
24
};
25
26
InfoMeasures compute_info_measures(const vector<double>& p_in, const vector<double>& q_in) {
27
    if (p_in.size() != q_in.size()) throw runtime_error("Size mismatch between p and q.");
28
    vector<double> p = normalize(p_in);
29
    vector<double> q = normalize(q_in);
30
31
    double H_P = 0.0, H_PQ = 0.0;
32
    for (size_t i = 0; i < p.size(); ++i) {
33
        double pi = p[i];
34
        if (pi == 0.0) continue; // term is 0 in both sums
35
        H_P  += -pi * safe_log(pi);
36
        H_PQ += -pi * safe_log(q[i]); // if q[i]=0 -> large penalty
37
    }
38
    double KL = H_PQ - H_P;
39
    return {H_P, H_PQ, KL};
40
}
41
42
int main() {
43
    // Example: P is near-uniform over 4 outcomes, Q is biased
44
    vector<double> P = {0.24, 0.26, 0.25, 0.25};
45
    vector<double> Q = {0.70, 0.10, 0.10, 0.10};
46
47
    try {
48
        InfoMeasures res = compute_info_measures(P, Q);
49
        cout.setf(ios::fixed); cout << setprecision(6);
50
        cout << "H(P)       = " << res.entropy_P << " nats\n";
51
        cout << "H(P,Q)     = " << res.cross_entropy << " nats\n";
52
        cout << "D_KL(P||Q) = " << res.kl_PQ << " nats\n";
53
        // Verify relation numerically
54
        cout << "Check: H(P)+KL = " << (res.entropy_P + res.kl_PQ) << "\n";
55
    } catch (const exception& e) {
56
        cerr << "Error: " << e.what() << "\n";
57
    }
58
    return 0;
59
}
60

1
#include <bits/stdc++.h>
2
using namespace std;
3
4
// Compute stable log-sum-exp of a vector
5
static inline double log_sum_exp(const vector<double>& z) {
6
    double m = *max_element(z.begin(), z.end());
7
    double sum = 0.0;
8
    for (double v : z) sum += exp(v - m);
9
    return m + log(sum);
10
}
11
12
// Given logits (batch x K) and integer labels (batch), compute mean loss and gradient dL/dz
13
struct SCEOutput { double loss; vector<vector<double>> grad; };
14
15
SCEOutput softmax_cross_entropy_with_grad(const vector<vector<double>>& logits,
16
                                          const vector<int>& labels) {
17
    size_t n = logits.size();
18
    if (n == 0) throw runtime_error("Empty batch");
19
    size_t K = logits[0].size();
20
    if (labels.size() != n) throw runtime_error("Label size mismatch");
21
    for (const auto& row : logits) if (row.size() != K) throw runtime_error("Jagged logits");
22
23
    vector<vector<double>> grad(n, vector<double>(K, 0.0));
24
    double loss_sum = 0.0;
25
26
    for (size_t i = 0; i < n; ++i) {
27
        int yi = labels[i];
28
        if (yi < 0 || (size_t)yi >= K) throw runtime_error("Label out of range");
29
        const auto& z = logits[i];
30
        double lse = log_sum_exp(z); // log denominator
31
        // Loss: -z[yi] + logsumexp(z)
32
        loss_sum += (-z[yi] + lse);
33
        // Gradient: softmax(z) - one_hot(yi)
34
        for (size_t k = 0; k < K; ++k) {
35
            double pk = exp(z[k] - lse); // softmax prob
36
            grad[i][k] = pk - (k == (size_t)yi ? 1.0 : 0.0);
37
        }
38
    }
39
40
    SCEOutput out;
41
    out.loss = loss_sum / static_cast<double>(n);
42
    // Average gradient over batch (common reduction)
43
    for (size_t i = 0; i < n; ++i)
44
        for (size_t k = 0; k < K; ++k)
45
            grad[i][k] /= static_cast<double>(n);
46
    out.grad = move(grad);
47
    return out;
48
}
49
50
int main() {
51
    // Example: 3 samples, 4 classes
52
    vector<vector<double>> logits = {
53
        {2.0, 1.0, 0.1, -0.5}, // true class 0
54
        {0.3, 1.2, -0.7, -0.2}, // true class 1
55
        {-1.0, 0.0, 0.5, 1.0}   // true class 3
56
    };
57
    vector<int> labels = {0, 1, 3};
58
59
    auto res = softmax_cross_entropy_with_grad(logits, labels);
60
    cout.setf(ios::fixed); cout << setprecision(6);
61
    cout << "Mean loss (nats): " << res.loss << "\n";
62
63
    cout << "Gradient (first sample):\n";
64
    for (double g : res.grad[0]) cout << g << ' ';
65
    cout << "\n";
66
    return 0;
67
}
68

1
#include <bits/stdc++.h>
2
using namespace std;
3
4
// Stable logistic loss: log(1 + exp(z)) - y*z
5
static inline double logistic_loss(double z, int y) {
6
    // For stability: use conditional form
7
    if (z >= 0) {
8
        return log1p(exp(-z)) + (1 - y) * z; // equals log(1+e^z) - y*z
9
    } else {
10
        return log1p(exp(z)) - y * z;        // avoids overflow when z << 0
11
    }
12
}
13
14
// Sigmoid for predictions (may underflow for extreme z but fine for reporting)
15
static inline double sigmoid(double z) {
16
    if (z >= 0) {
17
        double ez = exp(-z);
18
        return 1.0 / (1.0 + ez);
19
    } else {
20
        double ez = exp(z);
21
        return ez / (1.0 + ez);
22
    }
23
}
24
25
// Compute loss and gradient over batch X (n x d), labels y in {0,1}
26
struct LossGrad { double loss; vector<double> grad_w; double grad_b; };
27
28
LossGrad loss_and_grad(const vector<vector<double>>& X, const vector<int>& y,
29
                       const vector<double>& w, double b, double l2 = 0.0) {
30
    size_t n = X.size();
31
    size_t d = w.size();
32
    if (y.size() != n) throw runtime_error("Label size mismatch");
33
    for (const auto& xi : X) if (xi.size() != d) throw runtime_error("Feature size mismatch");
34
35
    vector<double> grad_w(d, 0.0);
36
    double grad_b = 0.0;
37
    double loss_sum = 0.0;
38
39
    for (size_t i = 0; i < n; ++i) {
40
        double z = inner_product(X[i].begin(), X[i].end(), w.begin(), 0.0) + b;
41
        loss_sum += logistic_loss(z, y[i]);
42
        double p = sigmoid(z);
43
        double err = p - static_cast<double>(y[i]); // derivative of loss wrt z
44
        for (size_t j = 0; j < d; ++j) grad_w[j] += err * X[i][j];
45
        grad_b += err;
46
    }
47
48
    // Average and add L2 regularization (if l2 > 0)
49
    double inv_n = 1.0 / static_cast<double>(n);
50
    for (size_t j = 0; j < d; ++j) grad_w[j] = grad_w[j] * inv_n + l2 * w[j];
51
    grad_b *= inv_n; // often bias is not regularized
52
    double reg = 0.5 * l2 * inner_product(w.begin(), w.end(), w.begin(), 0.0);
53
    double mean_loss = loss_sum * inv_n + reg;
54
55
    return {mean_loss, grad_w, grad_b};
56
}
57
58
int main() {
59
    // Create a simple linearly separable dataset in 2D
60
    vector<vector<double>> X; vector<int> y;
61
    std::mt19937 rng(42);
62
    normal_distribution<double> N1(-2.0, 1.0), N2(2.0, 1.0);
63
    for (int i = 0; i < 100; ++i) {
64
        double x0 = N1(rng), x1 = N1(rng); X.push_back({x0, x1}); y.push_back(0);
65
        double x2 = N2(rng), x3 = N2(rng); X.push_back({x2, x3}); y.push_back(1);
66
    }
67
68
    size_t d = 2;
69
    vector<double> w(d, 0.0);
70
    double b = 0.0;
71
    double lr = 0.1, l2 = 0.0;
72
73
    for (int it = 0; it < 200; ++it) {
74
        auto lg = loss_and_grad(X, y, w, b, l2);
75
        // Gradient descent update
76
        for (size_t j = 0; j < d; ++j) w[j] -= lr * lg.grad_w[j];
77
        b -= lr * lg.grad_b;
78
        if ((it+1) % 50 == 0) {
79
            cout.setf(ios::fixed); cout << setprecision(6);
80
            cout << "Iter " << (it+1) << ": loss = " << lg.loss << "\n";
81
        }
82
    }
83
84
    // Evaluate accuracy
85
    int correct = 0; int total = (int)X.size();
86
    for (int i = 0; i < total; ++i) {
87
        double z = inner_product(X[i].begin(), X[i].end(), w.begin(), 0.0) + b;
88
        int pred = sigmoid(z) >= 0.5 ? 1 : 0;
89
        correct += (pred == y[i]);
90
    }
91
    cout << "Training accuracy: " << (100.0 * correct / total) << "%\n";
92
93
    return 0;
94
}
95