Skip to content

English · Español

02 — Strides, vistas y broadcasting

Pruébalo — strides e índice plano

La estructura interna del ndarray de NumPy: un puntero a un buffer plano, una forma (shape), un vector de strides, un dtype, y unas flags. Esa estructura explica por qué arr.T es O(1) y np.ascontiguousarray(arr.T) es O(n). El broadcasting es un algoritmo de alineamiento de shapes que también necesitas dominar.


El ndarray, por dentro

Un ndarray de NumPy son cinco cosas:

ndarray = (
    data      : puntero a un buffer plano en memoria,
    shape     : tupla de ints (d_0, d_1, ..., d_{k-1}),
    strides   : tupla de ints en BYTES (s_0, s_1, ..., s_{k-1}),
    dtype     : tipo de elemento (fp32, int64, ...),
    flags     : metadatos (OWNDATA, C_CONTIGUOUS, F_CONTIGUOUS, WRITEABLE, ...),
)

El elemento (i_0, i_1, ..., i_{k-1}) vive en el offset en bytes:

\[ \text{offset}(i_0, \ldots, i_{k-1}) = \sum_{j=0}^{k-1} i_j \cdot s_j \]

Y ya está. Esa única fórmula es todo el modelo de memoria. Todo lo de abajo — vistas, transpuesta, broadcasting, indexación avanzada — es un corolario.

Ejemplo trabajado

import numpy as np
a = np.arange(12, dtype=np.int32).reshape(3, 4)
  • a.shape = (3, 4) — tres filas, cuatro columnas.
  • a.dtype.itemsize = 4 bytes (int32).
  • a.strides = (16, 4) — moverse una fila cuesta 16 bytes (4 elementos × 4 bytes), moverse una columna cuesta 4 bytes.
  • a.flags.C_CONTIGUOUS = True — row-major, el default.
  • a.flags.OWNDATA = Truea posee su buffer (él lo reservó).

El elemento a[2, 1] está en el offset 2*16 + 1*4 = 36 bytes dentro del buffer.

La transpuesta es gratis

b = a.T

¿Qué acaba de pasar? b comparte el mismo buffer (b.base is a es True), pero con shape y strides intercambiados:

  • b.shape = (4, 3).
  • b.strides = (4, 16).
  • b.flags.C_CONTIGUOUS = False (las filas ya no son contiguas en memoria).
  • b.flags.F_CONTIGUOUS = True (las columnas son ahora contiguas — orden Fortran).
  • b.flags.OWNDATA = False (b no posee el buffer; a sí).

El elemento b[1, 2] está en el offset 1*4 + 2*16 = 36 bytes — el mismo byte que a[2, 1]. La transpuesta es un re-etiquetado de los ejes; los bytes nunca se movieron.

Coste: O(1). Es solo una actualización de la struct.

Cuándo la transpuesta deja de ser gratis

En el momento en que pides una versión contigua:

c = np.ascontiguousarray(b)   # o b.copy() o np.asarray(b, order='C')

Ahora NumPy recorre b en su orden de strides no contiguo, copiando cada elemento a un buffer contiguo nuevo. Coste: O(n_elementos) en tiempo, O(n_elementos * itemsize) en memoria.

Por qué importa esto para código de IA: muchas rutinas BLAS / LAPACK (bajo np.linalg, np.matmul, np.dot) requieren entrada contigua. NumPy detecta el caso no contiguo y mete una copia oculta. El coste está oculto pero es real:

a = np.random.randn(1024, 1024).astype(np.float32)
b = np.random.randn(1024, 1024).astype(np.float32)

# Caso 1: ambos contiguos. El kernel de matmul corre directamente.
np.matmul(a, b)

# Caso 2: `a.T` es no contiguo. El kernel de matmul copia primero.
np.matmul(a.T, b)   # medible más lento

El lab 01 hace que Borja mida esto directamente. El número que verás en el i5-8250U para un float32 de 1024×1024 es aproximadamente: matmul contiguo ~30 ms, matmul no-contiguo-con-copia ~50 ms — la copia en sí es ~20 ms.

Vistas frente a copias: la tabla completa

Operación ¿Vista o copia?
a[1:3], a[::2], a[1:5:2] (slicing básico) Vista
a[1], a[1, 2] (indexación entera) Vista del subarray; escalar si todos los ejes están indexados
a[[1, 3, 5]] (indexación avanzada / por array) Copia
a[a > 0] (indexación booleana) Copia
a.T, a.transpose(...), np.swapaxes(a, ...) Vista
a.reshape(...) Vista si es posible (strides compatibles), si no copia
a.flatten() Copia (siempre)
a.ravel() Vista si es contiguo, copia en otro caso
np.ascontiguousarray(a) Copia si no es ya C-contiguo, no-op si lo es
a.copy() Copia (siempre)
a.astype(dtype) Copia si el dtype difiere; vista-o-copia si es el mismo dtype
np.asarray(a, dtype=X) No-op si a ya es un ndarray de dtype X, si no copia

Cómo comprobarlo en tiempo de ejecución:

  • arr.flags.OWNDATA — si es False, arr es una vista sobre el buffer de otro.
  • arr.base — el objeto original que la vista referencia, o None.
  • np.shares_memory(a, b) — definitivo (pero costoso — O(...) recorre los extents de strides de ambos arrays).

Trucos con strides (poder y peligro)

np.lib.stride_tricks.as_strided te deja construir un ndarray con shape y strides arbitrarios sobre un buffer existente. Así se implementan las ventanas deslizantes sin copiar:

from numpy.lib.stride_tricks import as_strided

a = np.arange(10, dtype=np.int32)
# Ventanas deslizantes de longitud 3:
window = as_strided(a, shape=(8, 3), strides=(4, 4))
# window[i, j] = a[i + j]

window es un array (8, 3) pero comparte el mismo buffer que a. Sin reserva de memoria.

El peligro: as_strided no hace comprobación de límites. Si especificas un extent shape × strides que pasa del buffer, lees memoria no inicializada. Si escribes en una vista con strides que se aliasa a sí misma, machacas datos. Trátala como solo-lectura.

El uso: NumPy moderno proporciona wrappers más seguros (np.lib.stride_tricks.sliding_window_view) que calculan los strides correctos automáticamente. Usa esos.

Broadcasting, formalizado

El broadcasting es el algoritmo de NumPy para operar sobre arrays de shapes diferentes. La regla es corta, las consecuencias son sutiles.

Las tres reglas

Dados dos shapes S_a = (a_0, a_1, ..., a_m) y S_b = (b_0, b_1, ..., b_n):

  1. Alinear a la derecha. Rellena el shape más corto con 1s a la izquierda hasta que ambos tengan el mismo número de ejes. Ej: (3,) vs (2, 4, 3)(1, 1, 3) vs (2, 4, 3).
  2. Compatibilidad por eje. Para cada eje, los dos tamaños de dim deben ser iguales o uno de ellos debe ser 1. Si ninguna, lanza ValueError.
  3. Shape resultado. Para cada eje, toma el máximo de los dos dims.

Ejemplos trabajados

(3,)         (2, 4, 3)    →    (2, 4, 3)
(N,)         (N, 1)       →    (N, N)   ← EL bug clásico
(B, 1, M)    (1, N, M)    →    (B, N, M)
(3, 4)       (3,)         →    ValueError  ← (3,) se hace (1, 3); (3, 4) necesita segundo dim 3 no 4. Espera, en realidad: alineados a la derecha, (3, 4) vs (3,) → (3, 4) vs (1, 3) → eje 1: 4 vs 3, NO compatible. ValueError.

El tercer ejemplo merece una pausa. np.array([1,2,3]) tiene shape (3,). Si quieres que broadcastee sobre las filas de un array (3, 4), tienes que darle shape (3, 1): np.array([1,2,3])[:, None]. El (3,) desnudo broadcastea sobre columnas.

El bug (N,) * (N,1)

y_pred = np.zeros((100,))      # shape (100,)
y_true = np.array([...]).reshape(100, 1)  # shape (100, 1)
err = y_pred - y_true          # shape (100, 100) !!

Broadcast alineado a la derecha: (100,)(1, 100), (100, 1) se queda. Ambos son ahora 2-D: (1, 100) vs (100, 1). Por eje: dim 0 es 1 vs 100 → broadcastea a 100; dim 1 es 100 vs 1 → broadcastea a 100. Resultado: (100, 100).

100×100 = 10.000 entradas. De estas, solo 100 son la diagonal que "habría sido correcta". Las otras 9.900 son términos cruzados. err.mean() promedia las 10.000 y devuelve el número equivocado.

Este bug es silencioso (sin error), incorrecto (devuelve un número que parece plausible) y ubicuo. El arreglo es siempre reshapear predicciones y objetivos al mismo shape antes de restar, siendo .reshape(-1, 1) o [:, None] el idioma estándar.

Por qué existe el broadcasting

No es una comodidad de Python; es una optimización de memoria. El broadcasting nunca materializa el array expandido. La expresión a + b donde a.shape = (N, 1) y b.shape = (1, N) actúa como si ambos se hubieran expandido a (N, N), pero solo reserva la salida (N, N). La "expansión" se hace con trucos de stride-cero por debajo.

np.broadcast_to(a, (N, N)) devuelve una vista con stride 0 sobre el eje broadcasteado — sin memoria reservada, pero el array se lee como si fuese (N, N).

Promoción de dtype

Cuando haces a + b y los dtypes difieren, NumPy promociona a un dtype común:

int32 + int64      → int64
int32 + float32    → float64    (sí, float64. la promoción de entero a float es generosa)
float32 + float64  → float64
int8 + bool        → int8

NumPy 2.0 cambió las reglas de promoción de escalares (NEP 50): np.float32(1) + 1.0 es ahora float32, no float64. La motivación: predictibilidad. Lee la NEP 50 una vez durante esta fase.

El patrón de bug que esto provoca: entrenas en fp32 a propósito, pero un numpy_array + python_float perdido se promociona a fp64 silenciosamente. Tu memoria se duplica. El lab no reproduce esto específicamente, pero la lista de escollos en PHASE_06_PLAN.md §5 lo señala.

Programación defensiva:

x = x.astype(np.float32, copy=False)   # castea o no-op
y = np.float32(0.5)                    # explícito
result = x + y                          # fp32 garantizado

copy=False es importante: astype(dtype, copy=True) siempre copia, incluso cuando el dtype ya coincide. copy=False es no-op cuando el dtype coincide.

Strides + broadcasting + dtype, combinados

Juntándolo todo: el coste de una expresión sobre ndarray depende de shape, strides, dtype y contigüidad además de la operación. Dos expresiones que parecen iguales pueden tener costes salvajemente distintos:

# (N, N) fp64 contiguo + (N,) fp64 → (N, N) fp64. Broadcast, sin copia. O(N²) trabajo, O(N²) escritura.
a + b

# (N, N) fp32 no contiguo (transpuesto) + (N,) fp32 → (N, N) fp32 salida contigua.
# Copia oculta de a antes de que corra el kernel. O(N²) trabajo + O(N²) copia.
a.T + b

El arreglo: perfilar (Fase 6 lab 03). La cura: predecir.

Recapitulación en un párrafo

Un ndarray de NumPy es (data, shape, strides, dtype, flags). El offset del elemento es Σ i_j · s_j — esa única fórmula explica por qué la transpuesta es O(1) (solo intercambiar shape y strides), por qué unas operaciones crean vistas y otras copian (depende de si un re-etiquetado solo-strides puede expresar el resultado) y por qué el broadcasting es rápido (magia de stride-cero, sin reserva del shape expandido). La regla de broadcasting es corta — alinear a la derecha, los dims deben coincidir o ser 1, el resultado es el máximo por pares — pero su modo de fallo silencioso ((N,) * (N,1) → (N,N)) es el bug de IA número 1. La promoción de dtype tiene sus propias sorpresas NEP-50. Domina esta página y una vasta clase de bugs futuros simplemente no te pasará.


Siguiente: 03-vectorization-and-profiling.md