Vectorization & Broadcasting — Tốc độ thật sự của NumPy

🎯 Mục tiêu

Sau bài này, bạn sẽ:

1. Hiểu tại sao vectorization nhanh hơn 10–100x so với vòng lặp Python thuần
2. Áp dụng broadcasting rules một cách tự tin — không còn đoán mò shape
3. Tư duy batch operations — nền tảng để viết code chạy được trên GPU

---

1. The Speed Problem — Python Loops vs NumPy

1.1 Bài toán: Chuẩn hóa 1 triệu giá trị

Giả sử bạn có 1 triệu số thực, cần chuẩn hóa min-max mỗi giá trị về đoạn [0, 1]:

$$x_{\text{norm}}=\frac{x-x_{min}}{x_{max}-x_{min}}$$

Có hai cách viết — và hiệu năng khác nhau tới 100 lần.

1.2 Cách 1: Vòng lặp Python thuần

import numpy as np
import time

# Tạo dữ liệu: 1 triệu giá trị ngẫu nhiên
data = np.random.rand(1_000_000)

# === CÁCH 1: Python for loop ===
start = time.perf_counter()

min_val = float('inf')
max_val = float('-inf')

# Tìm min, max bằng loop
for val in data:
    if val < min_val:
        min_val = val
    if val > max_val:
        max_val = val

# Chuẩn hóa bằng loop
result_loop = np.empty_like(data)
range_val = max_val - min_val
for i in range(len(data)):
    result_loop[i] = (data[i] - min_val) / range_val

loop_time = time.perf_counter() - start
print(f"Python loop: {loop_time:.4f}s")
# >>> Python loop: ~1.8 – 2.5s

1.3 Cách 2: NumPy vectorized

# === CÁCH 2: NumPy vectorized ===
start = time.perf_counter()

result_vec = (data - data.min()) / (data.max() - data.min())

vec_time = time.perf_counter() - start
print(f"NumPy vectorized: {vec_time:.6f}s")
# >>> NumPy vectorized: ~0.005 – 0.02s

# Kiểm tra kết quả giống nhau
print(f"Kết quả khớp: {np.allclose(result_loop, result_vec)}")
# >>> Kết quả khớp: True

print(f"Tốc độ gấp: {loop_time / vec_time:.0f}x")
# >>> Tốc độ gấp: ~100-200x

1.4 Tại sao chênh lệch khủng khiếp như vậy?

🔬 Giải thích kỹ thuật

Python loop chậm vì mỗi phần tử phải trả "thuế" interpreter:

Bước	Python Loop	NumPy Vectorized
Kiểm tra kiểu dữ liệu	Mỗi phần tử 1 lần	1 lần cho cả array
Gọi phép toán	1M lần gọi __sub__	1 lần gọi C function
Quản lý bộ nhớ	Tạo 1M Python objects	1 buffer liên tục
CPU cache	Cache miss liên tục	Sequential access = cache hit
SIMD instructions	❌ Không dùng được	✅ Xử lý 4-8 số cùng lúc

NumPy nhanh vì:

1. Compiled C loop — Vòng lặp chạy trong C, không qua interpreter
2. Contiguous memory — Dữ liệu nằm liền nhau → CPU cache hoạt động tối ưu
3. SIMD (Single Instruction, Multiple Data) — Một lệnh CPU xử lý 4–8 phần tử cùng lúc
4. Không tạo Python objects — Thao tác trực tiếp trên bộ nhớ raw

┌─────────────────────────────────────────────────────┐
│              PYTHON FOR LOOP                        │
│                                                     │
│  for i in range(1_000_000):                         │
│    ┌──────────┐  ┌──────────┐  ┌──────────┐        │
│    │ Type     │→ │ Dispatch │→ │ Execute  │→ next i │
│    │ Check    │  │ __sub__  │  │ float op │         │
│    └──────────┘  └──────────┘  └──────────┘        │
│    ↑ Lặp lại 1,000,000 lần qua interpreter ↑       │
└─────────────────────────────────────────────────────┘

┌─────────────────────────────────────────────────────┐
│              NUMPY VECTORIZED                       │
│                                                     │
│  data - data.min()                                  │
│    ┌──────────┐  ┌──────────────────────────────┐   │
│    │ 1x Type  │→ │ C loop + SIMD: xử lý 4-8    │   │
│    │ Check    │  │ phần tử mỗi chu kỳ CPU       │   │
│    └──────────┘  └──────────────────────────────┘   │
│    ↑ Overhead chỉ 1 lần, phần nặng chạy trong C ↑  │
└─────────────────────────────────────────────────────┘

---

2. Vectorization Mental Model

2.1 Tư duy "cả mảng", không phải "từng phần tử"

💡 Quy tắc vàng

Khi bạn thấy mình viết for i in range(len(array)) — DỪNG LẠI. Hỏi bản thân: "Phép toán này có thể diễn đạt bằng thao tác trên cả mảng không?"

Câu trả lời gần như luôn là CÓ.

Thay vì nghĩ:

"Với mỗi phần tử, trừ đi giá trị trung bình, rồi chia cho độ lệch chuẩn"

Hãy nghĩ:

"Trừ mảng cho mean, chia mảng cho std"

# ❌ Tư duy procedural (từng bước, từng phần tử)
for i in range(len(features)):
    features[i] = (features[i] - mean) / std

# ✅ Tư duy declarative (khai báo ý định trên cả mảng)
features = (features - mean) / std

2.2 Các pattern vectorization thường gặp

Pattern 1: Element-wise Operations

# Áp dụng cùng một phép toán lên mọi phần tử
squared = data ** 2                   # Bình phương tất cả
clipped = np.clip(data, 0, 1)        # Giới hạn về [0, 1]
activated = np.maximum(0, data)       # ReLU activation

Pattern 2: Reductions (gấp mảng thành 1 giá trị)

total = data.sum()                    # Tổng tất cả phần tử
avg = data.mean()                     # Trung bình
std = data.std()                      # Độ lệch chuẩn
idx_max = data.argmax()               # Vị trí phần tử lớn nhất

Pattern 3: Boolean Masking (lọc có điều kiện)

# Thay vì loop + if:
mask = data > 0.5                     # Mảng True/False
filtered = data[mask]                 # Chỉ lấy phần tử > 0.5
data[mask] = 0                        # Gán 0 cho phần tử > 0.5
count = (data > 0.5).sum()            # Đếm phần tử > 0.5

Pattern 4: Fancy Indexing

indices = np.array([0, 3, 7, 12])
selected = data[indices]              # Lấy phần tử tại các vị trí
data[indices] = -1                    # Gán giá trị tại các vị trí

2.3 Bảng chuyển đổi: Loop → Vectorized

Mục đích	❌ Python Loop	✅ Vectorized
Cộng 2 mảng	for i: c[i] = a[i] + b[i]	c = a + b
Tính tổng	s = 0; for x in a: s += x	s = a.sum()
Lọc điều kiện	for x in a: if x > 0: ...	a[a > 0]
Tìm max	m = a[0]; for x in a: ...	m = a.max()
Chuẩn hóa	for i: a[i] = (a[i]-μ)/σ	a = (a - μ) / σ
Nhân ma trận	for i: for j: for k: ...	C = A @ B
Đếm phần tử	c = 0; for x: if ... c+=1	(a > threshold).sum()

2.4 Bẫy .apply() của Pandas

⚠️ Pandas .apply() KHÔNG phải vectorization

Nhiều người lầm tưởng df.apply(func) là vectorized. Sai.

# ❌ Đây vẫn là Python loop ngụy trang!
df['norm'] = df['value'].apply(lambda x: (x - mean) / std)

# ✅ Đây mới là vectorized thật sự
df['norm'] = (df['value'] - mean) / std

Quy tắc: Nếu bạn truyền một Python function vào .apply(), Pandas vẫn gọi function đó từng dòng một — tốc độ tương đương for loop. Chỉ dùng .apply() khi logic quá phức tạp để vectorize.

---

3. Broadcasting Rules — Hướng dẫn đầy đủ

Broadcasting là cơ chế cho phép NumPy thực hiện phép toán giữa hai mảng có shape khác nhau mà không cần copy dữ liệu.

3.1 Ba quy tắc Broadcasting

📐 Ba quy tắc vàng

Quy tắc 1: Nếu hai mảng có số chiều (ndim) khác nhau, shape của mảng ít chiều hơn được thêm 1 vào đầu (bên trái) cho đến khi bằng nhau.

Quy tắc 2: Mảng có kích thước 1 ở một chiều nào đó sẽ hoạt động như thể được copy dọc theo chiều đó để khớp với mảng kia.

Quy tắc 3: Nếu ở bất kỳ chiều nào, hai mảng có kích thước khác nhau và không chiều nào bằng 1 → ValueError.

3.2 Quy tắc 1: Prepend 1s

Khi ndim khác nhau, NumPy tự động thêm 1 vào bên trái shape:

Ví dụ: mảng A shape (3, 4) + mảng B shape (4,)

Bước 1 – Cân bằng ndim:
   A: (3, 4)     ← 2 chiều, giữ nguyên
   B:    (4,)     ← 1 chiều, thêm 1 bên trái

   A: (3, 4)
   B: (1, 4)     ← Sau khi prepend

Bước 2 – Broadcasting:
   A: (3, 4)
   B: (3, 4)     ← dim 0: 1 → 3 (copy 3 lần)

Kết quả: (3, 4) ✓

import numpy as np

A = np.ones((3, 4))       # shape (3, 4)
B = np.array([1, 2, 3, 4])  # shape (4,)

C = A + B                  # shape (3, 4)
print(C)
# [[2. 3. 4. 5.]
#  [2. 3. 4. 5.]
#  [2. 3. 4. 5.]]

3.3 Quy tắc 2: Stretch dimensions of size 1

Chiều có kích thước 1 được "kéo dãn" (broadcast) để khớp:

Ví dụ: mảng A shape (3, 1) + mảng B shape (1, 4)

Đã cùng ndim = 2, so sánh từng chiều:
   dim 0: A=3, B=1 → B stretch thành 3
   dim 1: A=1, B=4 → A stretch thành 4

   A: (3, 1)  →  (3, 4)    ← cột được copy 4 lần
   B: (1, 4)  →  (3, 4)    ← hàng được copy 3 lần

Kết quả: (3, 4) ✓

Minh họa trực quan:

A (3,1):        B (1,4):          A + B (3,4):
┌───┐           ┌───┬───┬───┬───┐
│ 1 │           │ 10│ 20│ 30│ 40│
├───┤           └───┴───┴───┴───┘
│ 2 │     +                        =
├───┤
│ 3 │
└───┘

Broadcast A:    Broadcast B:       Kết quả:
┌───┬───┬───┬───┐  ┌───┬───┬───┬───┐  ┌───┬───┬───┬───┐
│ 1 │ 1 │ 1 │ 1 │  │ 10│ 20│ 30│ 40│  │ 11│ 21│ 31│ 41│
├───┼───┼───┼───┤  ├───┼───┼───┼───┤  ├───┼───┼───┼───┤
│ 2 │ 2 │ 2 │ 2 │  │ 10│ 20│ 30│ 40│  │ 12│ 22│ 32│ 42│
├───┼───┼───┼───┤  ├───┼───┼───┼───┤  ├───┼───┼───┼───┤
│ 3 │ 3 │ 3 │ 3 │  │ 10│ 20│ 30│ 40│  │ 13│ 23│ 33│ 43│
└───┴───┴───┴───┘  └───┴───┴───┴───┘  └───┴───┴───┴───┘

A = np.array([[1], [2], [3]])   # shape (3, 1)
B = np.array([[10, 20, 30, 40]])  # shape (1, 4)

C = A + B  # shape (3, 4)
print(C)
# [[11 21 31 41]
#  [12 22 32 42]
#  [13 23 33 43]]

3.4 Quy tắc 3: Incompatible shapes → Error

Ví dụ THẤT BẠI: mảng A shape (3, 4) + mảng B shape (3,)

Bước 1 – Cân bằng ndim:
   A: (3, 4)
   B: (1, 3)     ← prepend 1

Bước 2 – So sánh từng chiều:
   dim 0: A=3, B=1  → OK, B stretch thành 3
   dim 1: A=4, B=3  → ❌ FAIL! 4 ≠ 3 và không bên nào = 1

→ ValueError: operands could not be broadcast together
   with shapes (3,4) (3,)

A = np.ones((3, 4))
B = np.array([1, 2, 3])  # shape (3,) → (1, 3)

try:
    C = A + B
except ValueError as e:
    print(f"❌ Lỗi: {e}")
# ❌ Lỗi: operands could not be broadcast together
#    with shapes (3,4) (3,)

3.5 Thêm ví dụ: 3D Broadcasting

Ví dụ: Image batch (32, 224, 224, 3) + bias (3,)

Bước 1 – Cân bằng ndim:
   Image: (32, 224, 224, 3)
   Bias:                (3,) → (1, 1, 1, 3)

Bước 2 – Broadcasting:
   dim 0: 32 vs 1  → 1 stretch thành 32
   dim 1: 224 vs 1 → 1 stretch thành 224
   dim 2: 224 vs 1 → 1 stretch thành 224
   dim 3: 3 vs 3   → khớp

Kết quả: (32, 224, 224, 3) ✓
→ Bias được cộng vào mỗi pixel của mỗi ảnh!

3.6 Bảng tóm tắt Broadcasting

Shape A          Shape B          Kết quả
─────────────────────────────────────────────
(5,)             (5,)             (5,)        ← Cùng shape
(5,)             (1,)             (5,)        ← B stretch
(3, 4)           (4,)             (3, 4)      ← Rule 1 + 2
(3, 1)           (1, 4)           (3, 4)      ← Cả hai stretch
(3, 4)           (3,)             ❌ Error     ← Rule 3
(2, 3, 4)        (3, 4)           (2, 3, 4)   ← Rule 1
(2, 1, 4)        (1, 3, 1)        (2, 3, 4)   ← Nhiều stretch
(4, 3)           (5, 3)           ❌ Error     ← 4 ≠ 5

---

4. Real-World Broadcasting Examples

4.1 Feature Normalization — Trừ column means

Bài toán rất phổ biến trong ML: chuẩn hóa Z-score cho mỗi feature (cột).

# Dataset: 1000 samples, 5 features
X = np.random.randn(1000, 5) * 10 + 50  # shape (1000, 5)

# Tính mean và std theo cột (axis=0)
col_means = X.mean(axis=0)   # shape (5,)
col_stds = X.std(axis=0)     # shape (5,)

# Broadcasting magic:
#   X:         (1000, 5)
#   col_means: (5,) → (1, 5) → (1000, 5)   ← broadcast!
X_normalized = (X - col_means) / col_stds

# Kiểm tra: mỗi cột giờ có mean ≈ 0, std ≈ 1
print(f"Means sau chuẩn hóa:  {X_normalized.mean(axis=0).round(6)}")
print(f"Stds sau chuẩn hóa:   {X_normalized.std(axis=0).round(4)}")
# Means sau chuẩn hóa:  [-0.  0. -0. -0.  0.]
# Stds sau chuẩn hóa:   [1. 1. 1. 1. 1.]

Minh họa broadcasting:

X (1000, 5):              col_means (5,):
┌────┬────┬────┬────┬────┐
│ x₀₀│ x₀₁│ x₀₂│ x₀₃│ x₀₄│    [μ₀, μ₁, μ₂, μ₃, μ₄]
├────┼────┼────┼────┼────┤          ↓ broadcast ↓
│ x₁₀│ x₁₁│ x₁₂│ x₁₃│ x₁₄│    [μ₀, μ₁, μ₂, μ₃, μ₄]
├────┼────┼────┼────┼────┤    [μ₀, μ₁, μ₂, μ₃, μ₄]
│ ...│ ...│ ...│ ...│ ...│    ...
├────┼────┼────┼────┼────┤    [μ₀, μ₁, μ₂, μ₃, μ₄]
│xₙ₀ │xₙ₁ │xₙ₂ │xₙ₃ │xₙ₄ │
└────┴────┴────┴────┴────┘
Mỗi cột trừ đi mean của chính cột đó!

4.2 Neural Network — Cộng bias vector

Trong neural network, output của một layer là Y = X @ W + b:

batch_size = 32
input_dim = 128
output_dim = 64

X = np.random.randn(batch_size, input_dim)   # (32, 128)
W = np.random.randn(input_dim, output_dim)   # (128, 64)
b = np.random.randn(output_dim)              # (64,)

# Forward pass:
#   X @ W  → (32, 64)
#   b      → (64,) → (1, 64) → (32, 64)  ← broadcast!
Y = X @ W + b  # shape (32, 64) — mỗi sample cộng cùng bias

print(f"Output shape: {Y.shape}")
# Output shape: (32, 64)

💡 Tại sao bias broadcast?

Bias b có shape (64,) — mỗi neuron có 1 giá trị bias. Broadcasting tự động cộng cùng bias vào tất cả 32 samples trong batch. Không cần loop qua từng sample!

4.3 Discount Rates — Ma trận giao dịch

Áp dụng tỷ lệ chiết khấu khác nhau cho từng loại sản phẩm:

# 100 giao dịch × 5 loại sản phẩm
transactions = np.random.uniform(10, 500, size=(100, 5))  # (100, 5)

# Tỷ lệ chiết khấu cho mỗi loại sản phẩm
discounts = np.array([0.05, 0.10, 0.15, 0.08, 0.20])  # (5,)

# Broadcasting: (100, 5) * (5,) → (100, 5) * (1, 5) → (100, 5)
discount_amounts = transactions * discounts
final_prices = transactions - discount_amounts

# Hoặc gọn hơn:
final_prices = transactions * (1 - discounts)

print(f"Giao dịch gốc (sample 0):   {transactions[0].round(2)}")
print(f"Sau chiết khấu (sample 0):  {final_prices[0].round(2)}")

4.4 Distance Matrix — Khoảng cách giữa tất cả cặp điểm

# 5 điểm trong không gian 2D
points = np.random.randn(5, 2)  # shape (5, 2)

# Tính khoảng cách Euclidean giữa TẤT CẢ cặp điểm
# Trick: dùng broadcasting để tạo ma trận hiệu
#   points[:, np.newaxis, :] → (5, 1, 2)
#   points[np.newaxis, :, :] → (1, 5, 2)
#   Trừ nhau                 → (5, 5, 2)  ← broadcast!

diff = points[:, np.newaxis, :] - points[np.newaxis, :, :]
distances = np.sqrt((diff ** 2).sum(axis=-1))  # (5, 5)

print(f"Distance matrix shape: {distances.shape}")
print(f"Khoảng cách điểm 0 ↔ điểm 1: {distances[0, 1]:.4f}")

---

5. 🔥 GPU Connection — Tại sao vectorization là cánh cửa đến GPU

🔥 Tư duy GPU-friendly

Vectorized operations map trực tiếp sang GPU CUDA kernels. Mỗi phần tử trong mảng = một GPU thread. Broadcasting = GPU compiler tự động mở rộng dữ liệu để lấp đầy parallel grid.

Khi bạn viết Python loop cho array operations, GPU đắt tiền của bạn ngồi không trong khi Python interpreter chậm chạp xử lý từng phần tử một.

# ❌ GPU KHÔNG THỂ tăng tốc code này
for i in range(len(data)):
    result[i] = data[i] * weights[i] + bias

# ✅ GPU chạy TẤT CẢ phần tử ĐỒNG THỜI
result = data * weights + bias
# → Trên GPU: 1 triệu threads chạy cùng lúc!

Quy tắc: Nếu code NumPy của bạn không có for loop → nó có thể chạy trên GPU (với CuPy, JAX, hoặc PyTorch) gần như không cần sửa code.

---

6. 🧠 Misconception phổ biến

🧠 "Vectorization chỉ là về tốc độ" — SAI!

Vectorization còn mang lại CORRECTNESS (tính đúng đắn):

1. Không có off-by-one errors:

# ❌ Loop: dễ sai chỉ số
for i in range(len(data) - 1):  # Quên -1? Quên +1?
    diff[i] = data[i+1] - data[i]

# ✅ Vectorized: không có chỉ số nào để sai
diff = data[1:] - data[:-1]

2. Không có loop variable mutation bugs:

# ❌ Loop: biến total có thể bị ghi đè, quên reset
total = 0
for batch in batches:
    for x in batch:
        total += x  # Quên reset total giữa các batch?

# ✅ Vectorized: mỗi phép toán độc lập
totals = np.array([batch.sum() for batch in batches])

3. Dễ debug hơn vì operations là declarative:

# ❌ "Cho mỗi feature, trừ mean, chia std" — 5 dòng code, 3 biến tạm
# ✅ "Normalize all features" — 1 dòng, ý nghĩa rõ ràng
X_norm = (X - X.mean(axis=0)) / X.std(axis=0)

Vectorized code nói nó làm gì. Loop code che giấu ý định đằng sau mechanics.

---

7. ⚡ Fast Exercise — Column Normalization

⚡ Bài tập: Chuẩn hóa Z-score không dùng loop

Đề bài: Cho ma trận X có shape (100, 5) — 100 samples, 5 features. Chuẩn hóa mỗi cột về zero mean và unit variance (Z-score normalization).

Yêu cầu: KHÔNG dùng bất kỳ for loop nào.

Starter code:

import numpy as np

np.random.seed(42)
X = np.random.randn(100, 5) * np.array([10, 0.5, 100, 3, 25]) \
    + np.array([50, 2, 300, -5, 0])

# Feature scales rất khác nhau:
print("Means trước:", X.mean(axis=0).round(2))
print("Stds trước: ", X.std(axis=0).round(2))

# === YOUR CODE HERE ===
# X_normalized = ???
# ======================

# Kiểm tra:
# assert X_normalized.shape == (100, 5)
# assert np.allclose(X_normalized.mean(axis=0), 0, atol=1e-10)
# assert np.allclose(X_normalized.std(axis=0), 1, atol=1e-10)
# print("✅ Tất cả kiểm tra passed!")

---

Lời giải:

# Cách 1: Trực tiếp
X_normalized = (X - X.mean(axis=0)) / X.std(axis=0)

# Cách 2: Gán biến trung gian (dễ đọc hơn)
means = X.mean(axis=0)   # (5,) — mean mỗi cột
stds = X.std(axis=0)     # (5,) — std mỗi cột
X_normalized = (X - means) / stds

# Broadcasting:
#   X:     (100, 5)
#   means: (5,) → (1, 5) → (100, 5)
#   stds:  (5,) → (1, 5) → (100, 5)

assert X_normalized.shape == (100, 5)
assert np.allclose(X_normalized.mean(axis=0), 0, atol=1e-10)
assert np.allclose(X_normalized.std(axis=0), 1, atol=1e-10)
print("✅ Tất cả kiểm tra passed!")

---

8. 🪤 Gotcha — Broadcasting thầm lặng

🪤 Broadcasting "thành công" khi bạn muốn nó "thất bại"

Tình huống nguy hiểm nhất: Broadcasting tạo ra kết quả hợp lệ về mặt shape, nhưng sai về mặt logic.

a = np.array([[1], [2], [3]])   # shape (3, 1) — cột
b = np.array([10, 20, 30, 40])  # shape (4,)

# Bạn MUỐN cộng element-wise (3 phần tử + 3 phần tử)
# Nhưng b có 4 phần tử, không phải 3!

c = a + b  # Không báo lỗi! Kết quả shape (3, 4)
print(c)
# [[11 21 31 41]
#  [12 22 32 42]
#  [13 23 33 43]]

# Broadcasting vui vẻ tạo ma trận 3×4
# trong khi bạn muốn vector 3 phần tử!

Cách phòng tránh:

# 1. Luôn kiểm tra shape trước phép toán
assert a.shape == b.shape, f"Shape mismatch: {a.shape} vs {b.shape}"

# 2. Dùng np.add với casting kiểm soát chặt hơn
# 3. Print shape ở mỗi bước khi debug
print(f"a.shape={a.shape}, b.shape={b.shape}")

---

9. 📊 Performance Note — Nhân ma trận

📊 np.dot() vs @ vs np.matmul()

Syntax	Ý nghĩa	Khuyên dùng?
A @ B	Matrix multiplication	✅ Dùng cái này — ngắn gọn, Pythonic
np.matmul(A, B)	Giống @	Khi cần truyền như function argument
np.dot(A, B)	Dot product (hơi khác ở >2D)	⚠️ Tránh cho matrix — dùng cho vector
A * B	Element-wise multiply	Khác hoàn toàn với @!

A = np.random.randn(3, 4)
B = np.random.randn(4, 5)

# Ba cách viết — cùng kết quả cho 2D:
C1 = A @ B              # ✅ Ưu tiên dùng
C2 = np.matmul(A, B)    # OK
C3 = np.dot(A, B)       # OK cho 2D, cẩn thận với >2D

assert np.allclose(C1, C2) and np.allclose(C2, C3)

Quy tắc ngón tay cái:

• Element-wise: dùng * (phép nhân từng phần tử)
• Matrix multiply: dùng @ (tích ma trận)
• Đừng bao giờ nhầm lẫn hai cái này!

---

10. 🚫 Production Anti-pattern

🚫 Anti-pattern: iterrows() và apply(axis=1) trên triệu dòng

Trong production ML pipelines, đây là "sát thủ" hiệu năng số 1:

import pandas as pd

df = pd.DataFrame(np.random.randn(1_000_000, 5),
                  columns=['f1', 'f2', 'f3', 'f4', 'f5'])

# ❌ CHẬM KINH KHỦNG — ~30-60 giây cho 1M dòng
def compute_feature_slow(row):
    return row['f1'] * 2 + row['f2'] ** 2 - row['f3']

df['new_feat'] = df.apply(compute_feature_slow, axis=1)

# ❌ CÒN CHẬM HƠN — iterrows() tạo Series mỗi dòng
for idx, row in df.iterrows():
    df.loc[idx, 'new_feat'] = row['f1'] * 2 + row['f2']**2 - row['f3']

Cách đúng:

# ✅ Vectorized Pandas — ~0.01 giây cho 1M dòng
df['new_feat'] = df['f1'] * 2 + df['f2'] ** 2 - df['f3']

# ✅ Hoặc dùng .values để về NumPy nếu cần tốc hơn nữa
vals = df[['f1', 'f2', 'f3']].values  # shape (1M, 3)
df['new_feat'] = vals[:, 0] * 2 + vals[:, 1] ** 2 - vals[:, 2]

Phương pháp	Thời gian (1M dòng)	Tốc độ tương đối
iterrows()	~120s	1x (baseline)
apply(axis=1)	~30s	4x
Vectorized Pandas	~0.02s	6000x
NumPy .values	~0.01s	12000x

---

11. 🔬 Spot The Bug

🔬 Spot The Bug — Tìm lỗi trong code

Code sau đây có bug liên quan đến broadcasting. Bạn tìm được không?

import numpy as np

# Dữ liệu: 50 students × 4 subjects
scores = np.random.randint(0, 100, size=(50, 4))

# Trọng số cho mỗi môn: Toán=40%, Lý=30%, Hóa=20%, Anh=10%
weights = np.array([0.4, 0.3, 0.2, 0.1])

# Tính điểm trung bình có trọng số cho mỗi học sinh
weighted_avg = (scores * weights).sum()  # Bug ở đây!

print(f"Điểm TB mỗi HS: {weighted_avg}")

---

Lỗi: .sum() không có axis parameter → tính tổng tất cả phần tử, ra 1 số duy nhất thay vì 50 điểm trung bình.

Sửa:

# ✅ Đúng: sum theo axis=1 (theo hàng = theo mỗi học sinh)
weighted_avg = (scores * weights).sum(axis=1)  # shape (50,)

print(f"Shape: {weighted_avg.shape}")     # (50,)
print(f"HS đầu tiên: {weighted_avg[0]:.2f}")

Bài học: Khi dùng reduction (.sum(), .mean(), .std()), LUÔN chỉ định axis trừ khi bạn thực sự muốn gấp toàn bộ mảng thành 1 số.

---

12. 🎬 Scenario — Real ML Pipeline

🎬 Scenario: Feature Engineering cho bài toán dự đoán giá nhà

Bạn là Data Scientist tại một startup bất động sản. Dữ liệu có 100,000 giao dịch với 8 features:

import numpy as np

np.random.seed(42)
n_samples = 100_000

# Giả lập dữ liệu
data = {
    'area_m2': np.random.uniform(30, 200, n_samples),
    'rooms': np.random.randint(1, 6, n_samples),
    'floor': np.random.randint(1, 30, n_samples),
    'year_built': np.random.randint(1990, 2024, n_samples),
    'dist_center_km': np.random.uniform(0.5, 30, n_samples),
    'dist_metro_km': np.random.uniform(0.1, 10, n_samples),
    'has_parking': np.random.choice([0, 1], n_samples),
    'price_usd': np.random.uniform(50000, 500000, n_samples),
}

# Stack thành ma trận features
X = np.column_stack([data[k] for k in list(data.keys())[:-1]])  # (100000, 7)
y = data['price_usd']  # (100000,)

# === FEATURE ENGINEERING — TẤT CẢ VECTORIZED ===

# 1. Giá trên mét vuông (target engineering)
price_per_m2 = y / X[:, 0]  # Broadcasting: (100000,) / (100000,)

# 2. Tuổi nhà tính từ 2024
house_age = 2024 - X[:, 3]  # Broadcasting: scalar - (100000,)

# 3. Z-score normalization cho mọi feature
X_means = X.mean(axis=0)  # (7,)
X_stds = X.std(axis=0)    # (7,)
X_norm = (X - X_means) / X_stds  # Broadcasting: (100000,7) - (7,)

# 4. Interaction features: area × rooms
area_rooms = X[:, 0] * X[:, 1]  # (100000,)

# 5. Log transform cho features lệch phải
X_log = np.log1p(X[:, [0, 4, 5]])  # (100000, 3) — area, dist_center, dist_metro

print(f"Features chuẩn hóa — means:  {X_norm.mean(axis=0).round(6)}")
print(f"Features chuẩn hóa — stds:   {X_norm.std(axis=0).round(4)}")
print(f"Price/m²: min={price_per_m2.min():.0f}, max={price_per_m2.max():.0f}")
print(f"House age: min={house_age.min()}, max={house_age.max()}")

Không có for loop nào — tất cả 100,000 samples được xử lý song song qua vectorization và broadcasting.

---

13. 🛝 Playground — Code chạy thử

🛝 Playground: Vectorization Speedup + Broadcasting Demo

import numpy as np
import time

print("=" * 60)
print("DEMO 1: Vectorization Speedup")
print("=" * 60)

sizes = [1_000, 10_000, 100_000, 1_000_000]

for n in sizes:
    data = np.random.randn(n)

    # Python loop
    start = time.perf_counter()
    result_loop = [x ** 2 + 2 * x + 1 for x in data]
    loop_time = time.perf_counter() - start

    # NumPy vectorized
    start = time.perf_counter()
    result_np = data ** 2 + 2 * data + 1
    np_time = time.perf_counter() - start

    speedup = loop_time / np_time if np_time > 0 else float('inf')
    print(f"n={n:>10,}: Loop={loop_time:.4f}s, "
          f"NumPy={np_time:.6f}s, Speedup={speedup:.0f}x")

print()
print("=" * 60)
print("DEMO 2: Broadcasting Examples")
print("=" * 60)

# Example 1: Vector + Scalar
a = np.array([1, 2, 3, 4, 5])
print(f"\n[1] Vector + Scalar:")
print(f"    {a} + 10 = {a + 10}")

# Example 2: Matrix + Row vector
M = np.ones((3, 4))
v = np.array([10, 20, 30, 40])
print(f"\n[2] Matrix (3,4) + Row vector (4,):")
print(f"    Kết quả shape: {(M + v).shape}")
print(M + v)

# Example 3: Column × Row = Outer product
col = np.array([[1], [2], [3]])       # (3, 1)
row = np.array([10, 20, 30, 40])      # (4,)
print(f"\n[3] Column (3,1) * Row (4,) = Outer product (3,4):")
print(col * row)

# Example 4: Feature normalization
print(f"\n[4] Feature Normalization:")
X = np.random.randn(5, 3) * [10, 0.5, 100] + [50, 2, 300]
print(f"    Trước: means = {X.mean(axis=0).round(1)}")
X_norm = (X - X.mean(axis=0)) / X.std(axis=0)
print(f"    Sau:   means = {X_norm.mean(axis=0).round(6)}")
print(f"    Sau:   stds  = {X_norm.std(axis=0).round(4)}")

# Example 5: Broadcasting failure
print(f"\n[5] Broadcasting FAILURE:")
try:
    bad_a = np.ones((3, 4))
    bad_b = np.ones((3,))
    bad_c = bad_a + bad_b
except ValueError as e:
    print(f"    ❌ {e}")

print()
print("=" * 60)
print("DEMO 3: Performance — iterrows vs vectorized")
print("=" * 60)

# Simulate with NumPy (no pandas needed)
n = 100_000
data_matrix = np.random.randn(n, 3)

# "Loop" style
start = time.perf_counter()
result_loop = np.empty(n)
for i in range(n):
    result_loop[i] = data_matrix[i, 0] * 2 + data_matrix[i, 1] ** 2
loop_time = time.perf_counter() - start

# Vectorized
start = time.perf_counter()
result_vec = data_matrix[:, 0] * 2 + data_matrix[:, 1] ** 2
vec_time = time.perf_counter() - start

print(f"Loop (n={n:,}):       {loop_time:.4f}s")
print(f"Vectorized (n={n:,}): {vec_time:.6f}s")
print(f"Speedup:              {loop_time / vec_time:.0f}x")
print(f"Kết quả khớp:        {np.allclose(result_loop, result_vec)}")

---

14. Tổng kết

Checklist kiến thức bài này

Khái niệm	Bạn đã nắm?
Tại sao NumPy nhanh hơn Python loop 10–100x	⬜
Tư duy "whole array" thay vì "each element"	⬜
5+ patterns vectorization (element-wise, reduction, masking...)	⬜
Ba quy tắc broadcasting	⬜
Đọc ASCII shape diagram và predict kết quả	⬜
Broadcasting gotcha — silent shape mismatch	⬜
@ vs * — matrix multiply vs element-wise	⬜
Tại sao .apply() và iterrows() chậm	⬜

Quy tắc ghi nhớ

🧠 "Thấy loop trên array → Dừng lại, nghĩ vectorize"
📐 "Shape không khớp → Align phải, prepend 1, stretch 1s"
🔥 "Code không loop → Code chạy được trên GPU"
🪤 "Broadcasting thành công ≠ Logic đúng → Kiểm tra shape"

---

15. Bước tiếp theo

📚 Tiếp tục hành trình

Bài tiếp theo: Matrix Operations Intuition — Từ phép nhân ma trận đến hiểu bản chất của Linear Layers trong Neural Networks.

Bạn sẽ thấy: mọi thứ trong bài này (vectorization, broadcasting) là nền tảng để hiểu tại sao neural networks hoạt động — vì mỗi layer chỉ là một phép nhân ma trận + broadcasting bias!