Calificación de datos: estacionariedad y picos

Por José M. Requena Plens · Publicado el 21 jul 2026, actualizado el 26 jul 2026

Referencias: Bendat y Piersol 2010 Wald y Wolfowitz 1940 Rice 1945

Todo promedio de esta documentación - una PSD, un Leq, un balance GUM - supone que el registro es estacionario: que el proceso que lo genera no derivó mientras se medía. Bendat y Piersol dedican la Sección 10.3 de Random Data a calificar los registros antes del análisis, y phonometry.metrology implementa su núcleo cuantitativo: tests de tendencia y de estacionariedad libres de distribución con las regiones de aceptación del propio libro, y las estadísticas de Rice - cruces por nivel, frecuencia aparente, tasas y alturas de pico - que resumen el aspecto de un registro gaussiano calificado y delatan al que no lo es.

1. El test de inversiones de orden

Dada una secuencia de $N$ observaciones $x_{1}, \dots, x_{N}$ - estimaciones de parámetros, niveles por segmento, lo que sea - cuéntense los pares $i < j$ con $x_{i} > x_{j}$ . Cada par es una inversión de orden (reverse arrangement), y para observaciones independientes de una misma variable aleatoria su total $A$ tiene (B&P Ecs. (4.54)-(4.55))

μ_{A} = \frac{N ( N - 1 )}{4}, σ_{A}^{2} = \frac{N ( 2 N + 5 ) ( N - 1 )}{72},

sin ninguna suposición sobre la distribución de las $x_{i}$ . Una tendencia monótona empuja $A$ a un extremo (0 si la secuencia sube, $N (N - 1) /2$ si baja), así que la hipótesis de ausencia de tendencia se acepta al nivel de significación $α$ cuando $A$ cae dentro de una región bilateral - la Tabla A.6 de B&P, cuyas filas de $α = 0, 05$ reproduce exactamente trend_test y fija la suite de conformidad. El Ejemplo 4.4 del libro (veinte observaciones, $A = 86$ , aceptado entre 64 y 125) corre tal cual:

from phonometry import trend_test

values = [5.2, 6.2, 3.7, 6.4, 3.9, 4.0, 3.9, 5.3, 4.0, 4.6,
          5.9, 6.5, 4.3, 5.7, 3.1, 5.6, 5.2, 3.9, 6.2, 5.0]

res = trend_test(values)                  # B&P Ejemplo 4.4
print(res.statistic, res.bounds)          # 86, (64, 125)
print(res.trend_free, round(res.p_value, 3))   # True, 0.586
res.plot(language="es")

Dos secuencias de veinte observaciones dibujadas frente al índice de muestra: el Ejemplo 4.4 de Bendat y Piersol, que fluctúa en torno a cinco con ochenta y seis inversiones de orden y se acepta como sin tendencia, y las mismas fluctuaciones con una deriva ascendente añadida que sube de unos cinco a nueve, cuyas treinta y ocho inversiones de orden caen por debajo del límite inferior de aceptación de sesenta y cuatro y se rechazan como con tendencia al nivel del cinco por ciento

Mostrar el código de esta figura

import matplotlib.pyplot as plt
import numpy as np
from phonometry import trend_test

example = np.array([5.2, 6.2, 3.7, 6.4, 3.9, 4.0, 3.9, 5.3, 4.0, 4.6,
                    5.9, 6.5, 4.3, 5.7, 3.1, 5.6, 5.2, 3.9, 6.2, 5.0])
drifting = example + np.linspace(0.0, 4.0, example.size)   # deriva ascendente

res_flat = trend_test(example)
res_drift = trend_test(drifting)

fig, ax = plt.subplots(figsize=(10, 6))
index = np.arange(1, example.size + 1)
ax.plot(index, res_flat.values, "o-",
        label=f"Ejemplo 4.4: A = {res_flat.statistic}, aceptado")
ax.plot(index, res_drift.values, "s-",
        label=f"Deriva ascendente: A = {res_drift.statistic}, rechazado")
ax.set_xlabel("Índice de muestra")
ax.set_ylabel("Valor de la secuencia")
ax.legend()
plt.show()

El p-valor sale de la distribución nula exacta de $A$ (los conteos de inversiones de una permutación aleatoria), calculada hasta $N = 100$ - el alcance de la Tabla A.6 - y de la aproximación normal del libro más allá; el veredicto sigue la región tabulada del libro. .plot() dibuja la secuencia analizada frente a su índice de muestra con el conteo, la región de aceptación y el veredicto en la leyenda (con method="runs" marca además la mediana de la secuencia que clasifica cada valor).

2. Estacionariedad de un registro

El procedimiento de la Sec. 10.3.1.1 de B&P convierte el test de tendencia en un test de estacionariedad para una única historia temporal: divídase el registro en $N$ intervalos iguales lo bastante largos para ser independientes, calcúlese una media cuadrática por intervalo y pruébese esa secuencia. No hace falta saber nada del ancho de banda del registro, de su distribución ni de sus unidades, y el test funciona igual de bien sobre medias, valores eficaces o varianzas (statistic=). Una deriva de ganancia del 20 % - el escenario del Ejemplo 10.3 del libro - se detecta de inmediato, mientras que el mismo ruido sin deriva pasa:

import numpy as np
from phonometry import stationarity_test

fs = 8192.0
n = 1 << 16
noise = np.random.default_rng(42).standard_normal(n)

res = stationarity_test(noise, fs)        # 20 segmentos, medias cuadráticas
print(res.stationary, res.count, res.bounds)   # True, 91, (64, 125)

drifting = noise * np.linspace(1.0, 1.2, n)    # rampa de ganancia del +20 %
res = stationarity_test(drifting, fs)
print(res.stationary, res.count)          # False, 7  (tendencia al alza -> A bajo)
res.plot(language="es")

Veinte medias cuadráticas por segmento de dos registros de ruido: un registro estable cuyos valores fluctúan en torno a uno con un conteo de inversiones de orden de noventa y uno, aceptado como estacionario, y el mismo ruido con una rampa de ganancia del veinte por ciento cuyas medias cuadráticas por segmento suben de forma sostenida hasta uno coma cinco, dando solo siete inversiones de orden, rechazado como no estacionario al nivel del cinco por ciento

Mostrar el código de esta figura

import matplotlib.pyplot as plt
import numpy as np
from phonometry import stationarity_test

fs = 8192.0
n = 1 << 16
steady = np.random.default_rng(42).standard_normal(n)
ramp = np.random.default_rng(42).standard_normal(n) * np.linspace(1.0, 1.2, n)

res_steady = stationarity_test(steady, fs)
res_ramp = stationarity_test(ramp, fs)

fig, ax = plt.subplots(figsize=(10, 6))
index = np.arange(1, res_steady.n_segments + 1)
ax.plot(index, res_steady.segment_values, "o-",
        label=f"Ruido estable: A = {res_steady.count}, aceptado")
ax.plot(index, res_ramp.segment_values, "s-",
        label=f"Rampa de ganancia del +20 %: A = {res_ramp.count}, rechazado")
ax.set_xlabel("Índice de segmento")
ax.set_ylabel("Media cuadrática por segmento")
ax.legend()
plt.show()

El resultado registra la secuencia por segmento (segment_values, segment_times), el conteo, los bounds de la Tabla A.6, el p_value exacto y el veredicto stationary; .plot() dibuja la secuencia con el veredicto en la leyenda. El valor por defecto de 20 segmentos coincide con los ejemplos resueltos del libro - más segmentos resuelven derivas más rápidas, pero cada intervalo debe seguir siendo largo frente a las frecuencias más bajas del registro para que los valores sean independientes.

Merece la pena repetir dos salvedades del libro. Un registro puede ser no estacionario con media cuadrática estacionaria (un barrido de frecuencia, por ejemplo), así que pásese statistic="mean" o pruébense versiones filtradas por bandas cuando importe; y el test necesita que la tendencia sea lenta frente a la longitud del segmento, o se disuelve en la fluctuación aleatoria de los valores de segmento.

3. El test de rachas

El complemento clásico (tabulado en la tercera edición de Random Data; la distribución exacta es de Wald y Wolfowitz, 1940) clasifica cada valor como por encima o por debajo de la mediana de la secuencia y cuenta rachas de clasificación igual. Una tendencia o deriva lenta produce pocas rachas largas; la alternancia rápida produce demasiadas. method="runs" lo aplica con la distribución condicional exacta para cualquier $N$ - para veinte valores repartidos 10/10 la región de aceptación con $α = 0, 05$ es la clásica $(6, 15]$ :

import numpy as np
from phonometry import trend_test

rng = np.random.default_rng(3)
res = trend_test(rng.standard_normal(40), method="runs")
print(res.statistic, res.bounds, res.trend_free)
res.plot(language="es")   # la secuencia sobre su mediana con el veredicto

alternating = np.tile([1.0, -1.0], 10)   # 20 rachas: rechazado por el otro lado
print(trend_test(alternating, method="runs").trend_free)   # False

Dos paneles de secuencias clasificadas respecto a su mediana, puntos por encima en azul y por debajo en rojo. Izquierda: cuarenta observaciones gaussianas con veinte rachas dentro de la región de aceptación de catorce a veintisiete, veredicto sin tendencia. Derecha: una secuencia de veinte valores estrictamente alternante que también da veinte rachas, pero contra la región de aceptación de seis a quince para veinte valores, veredicto rechazado por alternancia demasiado rápida

El veredicto bilateral del test de rachas: 20 rachas no llaman la atención en 40 observaciones gaussianas (izquierda, aceptado), pero el mismo recuento en una secuencia alternante de 20 valores supera su cota superior (derecha, rechazado): demasiadas rachas es tan poco aleatorio como demasiado pocas.

Mostrar el código de esta figura

import matplotlib.pyplot as plt
import numpy as np
from phonometry import trend_test

rng = np.random.default_rng(3)
sequences = [rng.standard_normal(40), np.tile([1.0, -1.0], 10)]

# Una línea por secuencia: los valores analizados con el veredicto:
trend_test(sequences[0], method="runs").plot(language="es")
plt.show()

# Los dos veredictos lado a lado, clasificados respecto a la mediana:
fig, axes = plt.subplots(1, 2, figsize=(11, 4.5))
for ax, seq in zip(axes, sequences):
    res = trend_test(seq, method="runs")
    idx = np.arange(1, seq.size + 1)
    median = np.median(seq)
    above = seq > median
    ax.plot(idx, seq, "0.7", lw=0.7)
    ax.scatter(idx[above], seq[above], label="sobre la mediana")
    ax.scatter(idx[~above], seq[~above], color="r",
               label="bajo la mediana")
    ax.axhline(median, color="k", linestyle="--")
    lo, hi = res.bounds
    verdict = "sin tendencia" if res.trend_free else "rechazada"
    ax.set_title(f"r = {res.statistic} rachas, aceptación ({lo}, {hi}]: {verdict}")
    ax.set(xlabel="Índice de muestra", ylabel="Valor de la secuencia")
    ax.legend()
plt.show()

El test de inversiones de orden es el más potente de los dos frente a las tendencias monótonas que dominan la práctica (B&P Sec. 4.5.2), y por eso es el valor por defecto en todas partes; el de rachas añade sensibilidad al agrupamiento no monótono.

4. Cruces por nivel y frecuencia aparente

Para un registro gaussiano de media cero con autoespectro unilateral $G (f)$ , los resultados clásicos de Rice dan la tasa esperada de cruces por cero (ambas pendientes, B&P Ec. (5.195)) a partir de los momentos frecuenciales $m_{k} = \int f^{k} G (f) d f$ :

N_{0} = \frac{1}{π} \frac{σ _{v}}{σ _{x}} = 2 \frac{m _{2}}{m _{0}}, N_{a} = N_{0} e^{- a^{2} /2 σ_{x}^{2}},

donde $N_{a}$ es la tasa de cruces del nivel $a$ (Ec. (5.196)). $N_{0} /2$ es la frecuencia aparente del registro: un seno de 60 Hz cruza el cero 120 veces por segundo, el ruido paso bajo de ancho de banda $B$ da $N_{0} = 2 B / 3$ (una frecuencia aparente de $0, 58 B$ , Ejemplo 5.12) y una banda centrada en $f_{c}$ da $N_{0} = 2 f_{c}^{2} + B^{2} /12$ (Ejemplo 5.13). level_crossing_rate cuenta los cruces reales de cada nivel y pone al lado la curva de Rice, tomando los momentos del propio autoespectro de Welch del registro:

import numpy as np
from phonometry import level_crossing_rate

fs = 8192.0
t = np.arange(1 << 16) / fs
x = np.sin(2 * np.pi * 60.0 * t)          # un seno de 60 Hz...

res = level_crossing_rate(x, fs)
print(round(res.zero_crossing_rate, 1))   # ...tiene 120 ceros por segundo
print(round(res.apparent_frequency, 1))   # 60.0
res.plot(language="es")

Tasas de cruce por nivel medidas de un registro de ruido gaussiano limitado en banda entre ochocientos y mil doscientos hercios, dibujadas como puntos frente al nivel de cruce desde menos tres coma cinco hasta más tres coma cinco unidades de la señal en un eje de tasas logarítmico, cayendo desde unos dos mil cruces por segundo en el nivel cero hasta unos pocos por segundo a tres sigmas, con la curva exponencial de Rice pasando por todos los puntos medidos

Mostrar el código de esta figura

import matplotlib.pyplot as plt
import numpy as np
from phonometry import level_crossing_rate

fs = 20480.0
n = 1 << 19
rng = np.random.default_rng(0)
freqs = np.fft.rfftfreq(n, 1 / fs)        # ruido gaussiano limitado en banda:
spec = rng.standard_normal(freqs.size) + 1j * rng.standard_normal(freqs.size)
spec[(freqs < 800.0) | (freqs > 1200.0)] = 0.0
x = np.fft.irfft(spec, n)

res = level_crossing_rate(x, fs, levels=np.linspace(-3.5, 3.5, 29) * np.std(x))

fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(res.levels, res.rice_rates, label="Rice (Ec. 5.196)")
ax.plot(res.levels, res.rates, "o", label="Medida")
ax.set_yscale("log")
ax.set_xlabel("Nivel a [unidades de la señal]")
ax.set_ylabel("Cruces por segundo [1/s]")
ax.legend()
plt.show()

La curva de Rice vale para registros gaussianos, y ese es precisamente su segundo uso: unas tasas medidas que se apartan sistemáticamente de la curva son un cribado rápido de no gaussianidad (B&P Sec. 5.5.1.1) - el recorte aparece como cruces de nivel alto que faltan, la contaminación impulsiva como un exceso. Las tasas basadas en conteo necesitan el registro holgadamente sobremuestreado: un cruce entre dos muestras del mismo signo queda sin contar.

5. Picos: tasas, factor de irregularidad y alturas

Los mismos momentos fijan la tasa esperada de máximos locales $M = (1/2 π) (σ_{a} / σ_{v}) = m_{4} / m_{2}$ (Ec. (5.211)) y con ella el factor de irregularidad adimensional

r = \frac{N _{0}}{2 M} = \frac{m _{2}}{m _{0} m _{4}} \in (0, 1],

el número único que fija la distribución de las alturas de pico (B&P Sec. 5.5.4). Con $r = 1$ - datos de banda estrecha, un máximo por ciclo de cruce por cero - los picos siguen una distribución de Rayleigh: $Prob [pico > a] = e^{- a^{2} /2 σ_{x}^{2}}$ (Ec. (5.206)), la probabilidad de uno entre 3000 de un pico más allá de $4 σ$ del Ejemplo 5.14 de B&P. Cuando $r \to 0$ cada ciclo lleva cada vez más rizos encima, aparecen máximos negativos y las alturas de pico se acercan a la distribución de amplitud gaussiana simple. Entre medias, la mezcla de Rice (Ecs. (5.217)/(5.223)) interpola ambas, disponible como peak_exceedance() / peak_density() en el resultado:

import numpy as np
from phonometry import peak_statistics

fs = 8192.0
t = np.arange(1 << 16) / fs
x = np.sin(2 * np.pi * 60.0 * t)          # banda estrecha: r es en esencia 1

res = peak_statistics(x, fs)
print(round(res.irregularity_factor, 3))       # 0.997
print(res.peak_exceedance(4.0))                # exp(-8): aprox. 1 entre 3000
res.plot(language="es")   # excedencia empírica contra la mezcla de Rice (figura de abajo)

$Probabilidad de excedencia de la altura de pico de ruido gaussiano paso bajo en eje logarítmico frente a la altura de pico estandarizada desde menos dos coma cinco hasta cuatro coma cinco: la escalera empírica de medio millón de muestras sigue la curva de mezcla de Rice para un factor de irregularidad de cero coma siete cuatro seis, claramente por debajo del límite de Rayleigh discontinuo y por encima del límite gaussiano punteado, con una fracción visible de máximos negativos$ $Probabilidad de excedencia de la altura de pico de ruido gaussiano paso bajo en eje logarítmico frente a la altura de pico estandarizada desde menos dos coma cinco hasta cuatro coma cinco: la escalera empírica de medio millón de muestras sigue la curva de mezcla de Rice para un factor de irregularidad de cero coma siete cuatro seis, claramente por debajo del límite de Rayleigh discontinuo y por encima del límite gaussiano punteado, con una fracción visible de máximos negativos$

Mostrar el código de esta figura

import matplotlib.pyplot as plt
import numpy as np
from phonometry import peak_statistics
from phonometry.metrology.random_data import _rice_peak_exceedance

fs = 20480.0
n = 1 << 19
rng = np.random.default_rng(3)
freqs = np.fft.rfftfreq(n, 1 / fs)        # ruido paso bajo: r = sqrt(5)/3
spec = rng.standard_normal(freqs.size) + 1j * rng.standard_normal(freqs.size)
spec[freqs > 2000.0] = 0.0
x = np.fft.irfft(spec, n)

res = peak_statistics(x, fs)
peaks = res.peak_values
empirical = 1.0 - np.arange(1, peaks.size + 1) / peaks.size
z = np.linspace(-2.5, 4.5, 400)

fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(z, _rice_peak_exceedance(z, 1.0), "--", label="Rayleigh (r = 1)")
ax.plot(z, _rice_peak_exceedance(z, 0.0), ":", label="Gaussiana (r = 0)")
ax.plot(z, res.peak_exceedance(z),
        label=f"Rice (r = {res.irregularity_factor:.3f})")
ax.plot(peaks, empirical, drawstyle="steps-post", label="Empírica")
ax.set_yscale("log")
ax.set_ylim(1e-5, 1.5)
ax.set_xlabel(r"Altura de pico estandarizada $z = a/\sigma_x$")
ax.set_ylabel("Prob[pico > z]")
ax.legend()
plt.show()

Para una banda paso bajo ideal las formas cerradas dan $r = 5 /3 = 0, 745$ , y el valor medido aterriza encima. El factor de irregularidad es el puente estándar hacia la estimación de fatiga y de daño vibroacústico, donde selecciona la corrección de conteo de ciclos entre la hipótesis de Rayleigh de banda estrecha y las correcciones rainflow de banda ancha. Un aviso práctico: $m_{4}$ pondera el espectro por $f^{4}$ , así que el ruido de instrumentación de banda ancha muy por encima de la banda física infla $M$ y deprime $r$ sin avisar - limítese antes el registro a la banda físicamente significativa.

Dónde encaja esto

La calificación va antes de las estadísticas que calculan las demás páginas de esta sección: el intervalo de confianza chi-cuadrado de una PSD de Welch y las fórmulas de error aleatorio de los estimadores de correlación suponen que el registro es estacionario, igual que la idea misma de el Leq de una medición en Niveles. Cuando un registro no pasa el test, divídase en el cambio (la secuencia segment_values muestra dónde), analícense los trozos, o pásese a las vistas de tiempo corto - el espectrograma calibrado - que no suponen estacionariedad. Y cuando pasa, la maquinaria GUM puede propagar el error aleatorio restante de cada estimación promediada con la conciencia tranquila.

Qué cubre esta guía

Cubierto. La sección 4.5.2 con la Tabla A.6 de Bendat y Piersol (el test de tendencia por inversiones de orden, trend_test, su valor p exacto y por aproximación normal, y el Ejemplo 4.4), el test de rachas de Wald y Wolfowitz (1940) (trend_test(method="runs")), el procedimiento de estacionariedad por medias cuadráticas de segmento de la sección 10.3.1.1 con el Ejemplo 10.3 (stationarity_test), y las estadísticas de cruces por nivel y de picos de Rice (1945) de la sección 5.5 (level_crossing_rate, peak_statistics, el factor de irregularidad y sus límites Rayleigh y gaussiano para la altura de pico).

No cubierto. La sección 10.3 del libro se titula “calificación de datos” en un sentido más amplio y también cubre la clasificación, la validación y la edición de registros; phonometry.metrology implementa solo su núcleo cuantitativo, el test de estacionariedad por medias cuadráticas. Clasificar el tipo de un registro, validarlo frente a límites físicos o editar transitorios espurios siguen siendo pasos manuales tal como los describe el libro, fuera de este módulo.

Véase también

Referencia de la API: metrology.random_data.

Referencias

Bendat, J. S. y Piersol, A. G. (2010). Random data: Analysis and measurement procedures (4th ed.). Wiley. https://doi.org/10.1002/9781118032428Sección 4.5.2 con la Tabla A.6 (el test no paramétrico de tendencia por inversiones de orden, su media y varianza, los puntos porcentuales tabulados y el Ejemplo 4.4), Sección 10.3 (calificación de datos: clasificación, validación, edición; el procedimiento de estacionariedad por medias cuadráticas de segmento de 10.3.1.1 y el Ejemplo 10.3) y Sección 5.5 (cruces por nivel y valores de pico: Ecs. (5.195)-(5.196), (5.206), (5.211), (5.217)-(5.223) y Ejemplos 5.12-5.15). ISBN 978-0-470-24877-5. El test de rachas complementario aparecía en los puntos porcentuales de la distribución de rachas de la tercera edición.
Rice, S. O. (1945). Mathematical analysis of random noise. The Bell System Technical Journal, 24(1), 46-156. https://doi.org/10.1002/j.1538-7305.1945.tb00453.xLas derivaciones originales de las tasas esperadas de cruces por nivel y de máximos y de la distribución de picos del ruido gaussiano que presenta la Sección 5.5 de Bendat y Piersol (su Ref. 6; las Partes I-II están en el volumen 23, 1944).
Wald, A. y Wolfowitz, J. (1940). On a test whether two samples are from the same population. The Annals of Mathematical Statistics, 11(2), 147-162. https://doi.org/10.1214/aoms/1177731909La distribución condicional exacta del número de rachas, de la que se calculan las regiones de aceptación de rachas respecto de la mediana.

Creado y mantenido porJosé Manuel Requena Plens· GitHub· PyPI· Todos los proyectos