Capítulo 1. Conceptos de estadística y Modelamiento de datos

Estos apuntes para el curso AS750 se encuentran en http://www.das.uchile.cl/~mhamuy/courses/AS750/AS750.html
Sección 1. Incertezas en las mediciones

Lectura Obligatoria: Bevington (Capítulo 1)
Lectura Sugerida:

1.1 Definiciones básicas

• Error: Se define como la diferencia entre entre el valor observado o calculado y el valor verdadero.
  
• Exactitud: Medida de cuán cerca es la medición del valor verdadero.
  
• Precisión: Medida de cuán bien se determinó la medición, independientemente del valor verdadero.
• Errores sistemáticos: Errores que hacen que las mediciones se alejen del valor verdadero, independientemente de cuantas veces repitamos el experimento (ej. mala calibración del instrumento). Son una medida de la exactitud de la medición.
• Errores aleatorios: Errores debidos a fluctuaciones aleatorias que se disminuyen repitiendo el experimento. Son una medida de la precisión de la medición.
  
  Izquierda: mediciones precisas pero inexactas; Derecha: Mediciones exactas pero imprecisas

1.2 Distribución muestral y distribución madre

Si repetimos una medición muchas veces esperamos que los valores medidos se distribuyan alrededor del valor verdadero. Para el caso de un infinito número de mediciones obtendremos una distribución que determina la probabilidad de obtener un valor específico, la distribución madre. Esta distribución es hipotética ya que no podemos obtener infinitas medidas. En realidad uno obtiene una muestra finita, la distribución muestral.


Mediciones del largo de un bloque. La curva continua corresponde a una Gaussiana centrada en 19.9 cm y σ=0.52 cm. La curva segmentada corresponde a la distribución Gaussiana madre centrada en 20.0 cm y σ=0.50 cm.

El área bajo la curva segmentada entre x y x+dx da el número de eventos esperado en ese intervalo en una muestra de 100 mediciones. Esta área dividida por el área total es la probabilidad P(x)dx de obtener una medición entre x y x+dx. P(x) es la función de densidad de probabilidad.

1.3 Media, moda y mediana

• Valor medio de la muestra:
  
• Valor medio de la población madre:
  
• Mediana de la población madre:
  

  La mediana divide la distribución en dos áreas iguales.

• Moda de la población madre:
  

  La moda es el valor más probable.

  
  Distribución asimétrica de probabilidad que ilustra la media, la moda y la mediana de la variable x.

1.4 Desviaciones



El promedio de las desviaciones se anula debido a la definición de la media:



La desviación media se define por:

y corresponde a una medida de la dispersión de las mediciones alrededor de la media.

La varianza se define como el promedio de los cuadrados de las desviaciones:



La desviación estándar σ se define como la raíz cuadrada de la varianza.

La desviación estándar de la muestra, s, se define como:



Se debe usar el factor N-1 en vez de N ya que la media se determina a partir de los datos.

La media μ y la desviación estándar σ son los parámetros que necesitamos para describir una distribución simétrica. La media es la mejor estimación del valor verdadero. La desviación estándar es una medida de la incerteza de nuestras mediciones. Para describir distribuciones asimétricas se requieren parámetros adicionales.

1.5 Media y desviación estándar de distribuciones discretas

Si la función de probabilidad P(x) es conocida, el valor medio se calcula como:



en que N es el número de mediciones y n es el número de bins.

De igual modo, la desviación estándar es:



En general, el valor esperado para cualquier función f(x) es:



1.6 Media y desviación estándar de distribuciones continuas

Cuando la distribución es continua se reemplaza la sumatoria por una integral para calcular los valores esperados.





La distribución madre no se puede determinar. El valor medio <x> de la muestra es la mejor estimación de μ mientras que la desviación estándar s de la muestra es la mejor estimación de σ.

Sección 2. Distribuciones de probabilidad

Lectura Obligatoria: Bevington (Capítulo 2)
Lectura Sugerida:

2.1 La distribución binomial

Si p es la probabilidad de un cierto evento individual (ej: p=1/2 para que una moneda salga cara o p=1/6 para que un dado de 6) y n es el número de lanzamientos de una moneda o un dado, la probabilidad de obtener x de dichos eventos (x≤n) está dada por la distribución binomial:



La distribución binomial es una distribución discreta, definida sólo para números enteros positivos de x.

Valor medio:



Varianza:



Ejemplo.- Se lanzan 10 monedas al aire


Distribución binomial para μ=5 y p=0.5.

Ejemplo.- Se lanzan 10 dados


Distribución binomial para μ=10/6 y p=1/6.

2.2 La distribución de Poisson

La distribución de Poisson corresponde al caso particular de la distribución binomial cuando p<<1.



Esta distribución es la apropiada para experimentos de conteo de eventos (ej: fotones) por intervalo de tiempo.

Valor medio:



Varianza:



Notar que la desviación estándar σ es la raíz cuadrada de la media, por lo cual la distribución de Poisson está descrita por un solo parámetro.

Notar que, si bien el error absoluto al medir N cuentas aumenta como √N, el error relativo disminuye como 1/√N.

Ejemplo.- Medición del número de cuentas de un isótopo radioactivo


Histograma de detección de cuentas en intervalos de 2 segundos. Se grafica con línea continua la distribución de Poisson para μ=1.69.

Histograma de detección de cuentas en intervalos de 15 segundos. Se grafica con línea continua la distribución de Poisson para μ=11.48.

2.3 La distribución de Gauss (normal)

La distribución de Gauss corresponde a la distribución de Poisson cuando μ se hace grande.



Distribución de probabilidad de Gauss.

Esta distribución da la probabilidad de obtener el valor x a partir de una distribución madre con parámetros μ y σ, correspondientes a la media y desviación estándar. Esta es una función continua (definida para cualquier valor de x) y normalizada:



Esta distribución es simétrica alrededor de μ. El ancho de la curva está caracterizado por σ:



La integral de la funcion de probabilidad de Gauss esta tabulada en la tabla C.2. Hay un 68% y 95% de probabilidad de que una medición caiga dentro de 1σ y 2σ de la media.

Ancho total a mitad de altura:




2.4 La distribución de Lorentz

Esta distribución es apropiada para datos de comportamiento de renonancia.



Distribuciones de probabilidad Gaussiana y Lorentziana.

Esta función es simétrica alrededor de la media μ. Su ancho está dado por su ancho medio Γ. Esta función tiene alas más extendidas que las de la distribución de Gauss. Es una función normalizada:



La desviación estándar diverge:


Sección 3. Análisis de errores

Lectura Obligatoria: Bevington (Capítulo 3)
Lectura Sugerida:

3.1 Errores internos versus errores externos

Cuando la estimación de los errores utiliza suposiciones sobre la distribución madre (por ejemplo al asumir errores en errores del instrumento, tal como la incerteza al medir el largo de una mesa con una regla graduada) se obtiene una medida del error externo.

Una estimación de los errores a partir de la varianza es una medida empírica de la dispersión de las mediciones y del acuerdo relativo de los datos. Esta estimación a partir de los mismos datos se denomina error interno.

Si las suposiciones sobre los errores instrumentales son correctas, los errores internos y externos debieran coincidir.

3.2 Propagación de errores

A menudo uno necesita saber el error en una cantidad que es función de una o más cantidades medidas. Cuál es el error en la cantidad de interés?

Por ejemplo, si estuviésemos interesados en el volumen V de una caja de largo L, ancho W, y altura H, cuyos errores fueran ΔL=L-L₀, ΔW=W-W₀, ΔH=H-H₀, podriamos estimar el error en el volumen ΔV=V-V₀ expandiendo V en series de Taylor:



Evaluando las derivadas parciales:



podemos calcular ΔV a partir de los errores ΔL, ΔW, ΔH.

En el caso general de una cantidad x que depende de dos cantides medidas u, v:



la incerteza en x se puede estimar a partir de la dispersión en los valores individuales x_i:



La varianza en x es:



Las desviaciones x-x_i se pueden obtener a partir de las desviaciones en las cantidades medidas u-u_i, v-v_i:

con lo cual:



Los primeros dos términos se pueden expresar en términos de las varianzas de u y v:




El tercer término se puede reescribir introduciendo la covarianza entre las variables u y v:



con lo cual obtenemos la fórmula de propagación de errores:



Si las fluctuaciones en u y v no están correlacionadas esperamos que el término cruzado se cancele con lo cual:



Ejemplo de fluctuaciones con y sin correlación:


Panel superior: 1000 mediciones no correlacionadas del ancho y largo de un bloque. Panel inferior: el volumen (LxWxH) y el area del bloque (LxW) calculados a partir de mediciones no correlacionadas de L, W y H resultan altamente correlacionadas.

Ejemplo 1.-



Ejemplo 2.-



Ejemplo 3 (potencias).-



Ejemplo 4 (exponentes).-



Ejemplo 5 (logaritmos).-



Ejemplo 6 (magnitudes).-



Regla de oro: un error de 0.01 mag corresponde a un error de 1% en el flujo medido

Sección 4. Estimaciones de promedios y errores

Lectura Obligatoria: Bevington (Capítulo 4)
Lectura Sugerida:

4.1 Método de cuadrados mínimos

Supongamos que en un cierto experimento tenemos N observaciones aleatoriamente extraidas de la población madre (un muestra infinita distribuida de acuerdo a la distribución madre). Suponiendo que la distribución madre es una Gaussiana con media μ y desviación estándar σ, la probabilidad dQ_i de realizar una observacion x_i en el intervalo dx es:



En principio no conocemos la media μ y debemos determinarla a partir de los datos. Si denominamos esta cantidad como μ', la probabilidad está dada por:



La probabilidad de obtener el set de N observaciones es el producto de las probabilidades individuales:



o bien:



El método de la máxima probabilidad dice que lo más probable es que los datos provengan de una distribución con μ=μ', aquel valor de μ' que maximice P(μ'). Maximizar P(μ') es equivalente a maximizar el argumento X de la exponencial:



Para encontrar μ' debemos buscar el valor que anule esta derivada:




Si todas las mediciones de x se hicieron con la misma precisión, σ es constante:



El valor que maximiza la probabilidad corresponde a la media <x>.

Cuál es la incerteza σ asociada a la media μ'? Cada medición x_i tiene asociada una incerteza σ, de manera que la media debe tener un error menor a cada medición individual. Usando la fórmula de propagación de errores:



La derivada parcial es:



Como σ_i=σ=constante,



El error en la media μ disminuye con la raíz cuadrada del número de mediciones. La desviación estándar de la población madre se puede estimar internamente a partir de los mismos datos:



con lo cual podemos obtener la incerteza en la media μ:



Qué hacemos cuando las mediciones individuales x_i se realizan con distinta precisión σ_i? La probabilidad de obtener los N datos es:



Usando el método de la máxima probabilidad debemos maximizar esta probabilidad, lo cual equivale a maximizar el argumento de la exponencial:



El valor más probable es la media ponderada de los datos:



Para calcular el error en la media calculamos nuevamente la derivada parcial:



con lo cual:



4.2 Test de χ²

En esta sección describiremos el test de χ², el cual nos permite determinar si un cierto modelo es una buena descripción de los datos.

Si realizamos N mediciones de la cantidad x obtendremos N valores de x_i. Luego podemos contar cuantas veces se repite cada valor de x_j. Si denominamos h(x_j) el número de eventos en cada bin podemos construir un histograma de j=1→n bins, tal como se muestra en rojo en esta figura:


Histograma de 100 mediciones construido de una distribución Gaussiana con μ=5 y σ=1. La distribución madre se muestra en azul. Las curvas de puntos negros representan la distribución Poissoniana de eventos en cada bin, calculada en base a la distribución Gaussiana.

Si P(x_j) es la probabilidad de observar el valor x_j entonces debemos esperar y(x_j)=NP(x_j) eventos en dicho bin (curva azul).

Si repitiéramos este experimento muchas veces (k=1→n_k), obtendríamos muchas mediciones de h_k(x_j) para cada bin y podemos encontrar la distribución p_j(y_k) que describe la probabilidad de obtener el valor h_k(x_j). Como el experimento de medir h_k(x_j) consiste en medir eventos discretos, esperamos que p_j(y_k) corresponda a una distribución de Poisson (curvas punteadas negras). Según la estadística de Poisson la desviación estándar esperada en cada bin debe estar dada por σ_j(h)=√y(x_j).

En este ejemplo estamos asumiendo que la distribución madre es una Gaussiana pero en un experimento real no conocemos la forma de la distribución madre, por lo cual la desviación estándar distribución Poissoniana en cada bin se debe calcular a partir de los mismos datos, i.e., σ_j(h)=√NP(x_j). Esto se ilustra como curvas punteadas negras en esta figura:


Histograma de 100 mediciones construido de una distribución Gaussiana con μ=5 y σ=1. La distribución madre se muestra en azul. Las curvas de puntos negros representan la distribución Piossoniana de eventos en cada bin, calculada en base a la muestra de datos observada.

Usando las definiciones anteriores se define el parámetro χ² (chi cuadrado) como:



En un experimento real podemos estimar la incerteza en cada bin a partir de los datos:



con lo cual la definición de χ² se simplifica a:



El numerador es una medida de la dispersión en las observaciones mientras que el denominador es una medida de la dispersión esperada.

• Si el acuerdo entre los datos y el modelo fuera perfecto obtendríamos χ²=0.
• Si la dispersión de los datos correspondiera a la dispersión esperada, esperaríamos una contribución de una unidad por cada bin, χ²≅n para todo el histograma, o bien
• Sin embargo, en realidad obtendremos diferencias con el valor esperado y debemos saber cuál es la probabilidad de medir un cierto valor de χ² a partir de una muestra aleatoria de datos.

Se define el parámetro chi cuadrado reducido como χ²_ν=χ²/ν.

• El valor esperado de χ²_ν es 1.
• Valores de χ²_ν mucho mayores que 1 indican grandes diferencias entre los datos y lo esperado.
• Valores de χ²_ν muy pequeños tampoco son aceptables.

La tabla C.4 da la probabilidad de que una muestra aleatoria de datos extraída de una cierta distribución de probabilidad arroje un valor de χ² tan grande o mayor que el valor observado en un experimento con ν grados de libertad. Si la probabilidad es cercana a uno entonces la distribución de probabilidad asumida es un buen modelo. De este modo la estadística χ² nos permite juzgar si el modelo es una buena descripción de los datos.

En esta sección se asumió una distribución de Poisson para modelar la dispersión en cada bin. El test χ² se puede generalizar para comparar observaciones de cualquier experimento a las predicciones de cualquier modelo.

Sección 5. Ajuste de cuadrados mínimos a una línea recta

Lectura Obligatoria: Bevington (Capítulo 6)
Lectura Sugerida: Numerical Recipes (Secciones 15.0, 15.1, 15.2)

5.1 Método de cuadrados mínimos

En muchas experimentos se mide más de una variable, es decir en vez de medir sólo la cantidad x, realizamos N mediciones de (x_i,y_i). Nuestro objetivo es encontrar una función y=y(x) que describa adecuadamente la relación entre estas dos variables. En este capítulo asumiremos una simple relación lineal (el "modelo"):


Ejemplo 1.- Medición de la diferencia de potencial en función de la posición de un alambre conductor:


Diferencia de potencial como función de la posición de un alambre conductor. La recta es un ajuste de cuadrados mínimos a los datos.

En este ejemplo el estudiante estima que la incerteza en las distancias es despreciable, mientras que para el voltaje asume un error de 0.05 V.

Ejemplo 2.- Mediciones de cuentas de una fuente radioactiva en función de la distancia a la fuente:


Número de cuentas en intervalos constantes de tiempo de una fuente radioactiva como función del inverso del cuadrado de la distancia a la fuente. La recta es un ajuste de cuadrados mínimos a los datos.

En este ejemplo el estudiante estima que la incerteza en las distancias es despreciable, y para el número de cuentas asume errores variables dados por la estadística de Poisson.

Em ambos ejemplos queremos ajustar una relación lineal entre dos variables:



Nuestro objetivo es determinar los coeficientes a y b que mejor reproducen las mediciones. Para esto contamos con N mediciones (x_i,y_i). Se define la desviación Δy_i como la diferencia entre el valor medido y el valor calculado del modelo:



Queremos calcular los coeficientes a₀ y b₀



Asumiendo que cada medición y_i se rige por una distribución Gaussiana centrada en y₀(x_i) y desviación estándar σ_i, la probabilidad de realizar dicha medición es:



La probabilidad de obtener el set de N mediciones es:



Podemos calcular esta probabilidad para cualquier par de parámetros (a,b):



El método de máxima probabilidad sostiene que los valores más probables de (a,b) son aquellos que maximizan esta probabilidad. Esto es equivalente a maximizar le suma de la exponencial, la cual se define como el parámetro χ²:



Los valores de a y b que minimizan χ² son aquellos para los cuales las derivadas de χ² se anulan:




Estas restricciones constituyen dos ecuaciones lineales para las incógnitas a y b:




Usando el método de los discriminantes:




en que:



Para el caso especial en que todas las incertezas son iguales (σ=σ_i):





5.2 Incerteza en los parámetros

Para calcular la incerteza en los parámetros a y b usaremos la fórmula de propagación de errores:



Para esto necesitamos las derivadas parciales de a y b con respecto a cada y_j que contribuyó al cálculo de a y b:




Con esto obtenemos la varianza en a y b:








Para el caso especial en que σ=σ_i:


Sección 6. Ajuste de cuadrados mínimos a un polinomio

Lectura Obligatoria: Bevington (Capítulo 7)
Lectura Sugerida: Numerical Recipes (Sección 15.4)

6.1 Solución por el método de determinante

En esta sección usaremos una relación más compleja para modelar los datos, en particular una función polinomial del tipo:



que se puede escribir como:



En este modelo hay m coeficientes a₁, a₂, ... a_m a ajustar. Podemos generalizar el método incluyendo funciones arbitrarias f_k(x), siempre que éstas NO incluyan a los coeficientes a_k.



La probabilidad combinada de obtener los N pares de mediciones (x_i,y_i) es:



Podemos definir el parámetro χ² como la suma:



Usando nuevamente el método de máxima probabilidad, debemos encontrar los parámetros a_k que maximizan esta probabilidad:



Esta condición arroja m ecuaciones acopladas y lineales en los coeficientes a_k:



La solución a estas m ecuaciones se puede obtener por el método de determinante.

Ejemplo.- Medición del voltaje de una termocupla como función de la temperatura. Se asume que los errores en la temperatura son despreciables mientras que las incertezas en los voltajes se consideran constantes (0.05 mV).


Voltaje de una termocupla versus temperatura. La curva corresponde a un polinomio de segundo orden obtenido por el método de cuadrados mínimos.

La relación V/T es aproximadamente lineal pero una inspección más cercana de los datos revela una relación de la forma V = a₁ + a₂ T + a₃T². En este caso la función de χ² es:



La suma se realiza hasta i=21 que corresponde al número de mediciones. El ajuste de cuadrados mínimos arroja χ²=26.6, o bien χ²_ν=χ²/18=1.5. La tabla C.4 indica que la probabilidad de obtener χ²_ν igual o mayor a 1.5 es de 8.8%. Esta probabilidad indica que el ajuste es razonable.

6.2 Solución por el método de matrices

Para poder determinar los errores en los parámetros del ajuste utilizaremos en esta sección el método de las matrices que resulta elegante y conveniente. Recordemos que estamos tratando de ajustar la función:


La solución para los coeficientes a_k se obtiene de m ecuaciones lineales en los coeficientes:



Estas ecuaciones se pueden escribir de manera matricial:



en que los elementos de las matrices α y β son:




y los elementos de la matriz a son los coeficientes a_k del ajuste. Para encontrar la matriz a debemos multiplicar ambos lados de la ecuación por la matriz ε, la inversa de α:



Esto da:



lo que se puede expresar como:



Así, el problema de encontrar los coeficientes a_l se reduce a invertir la matriz α.

La covarianza σ_ajal entre los parámetros a_j y a_l es:



Las derivadas parciales son:



con lo cual la covarianza queda:



Se desprende entonces que los elementos de la matriz ε contiene las covarianzas entre los distintos coeficientes. La matriz ε se denomina la matriz de covarianza. Sus elementos permiten calcular las incertezas en los coeficientes ajustados y el grado de covarianza entre diferentes coeficientes.

Ejemplo.- Distribución angular de una fuente de rayos γ.


Distribución angular de una fuente de rayos γ. La curva corresponde a un ajuste de cuadrados mínimos a una serie de polinomios de Légendre hasta L=4.

En este ejemplo se usaron polinomios de Légendre para realizar el ajuste.

6.3 Linealización de funciones no lineales

En todos los ejemplos anteriores hemos asumido que la función de ajuste es lineal en los parámetros a_k. Esta característica lineal es lo que permite obtener una solución analítica a los coeficientes. Cuando la función pierde esta característica, en general no se pueden obtener soluciones analíticas y la solución debe buscarse numéricamente. Hay casos de funciones no lineales, sin embargo, que se pueden linealizar mediante un cambio de variables.

Por ejemplo, la función siguiente es claramente no lineal en el coeficiente b:



Mediante un cambio de variables se puede linealizar:



y la solución se puede obtener con los métodos anteriores. Hay que ser cuidadosos con los errores en la medición de la variable y. Se debe usar la fórmula de propagación de errores para calcular las incertezas en la nueva variable ln(y).

La linealización de la función original se puede ver gráficamente:


Gráfico de la función y=1000×e^-5x. Las barras de error están dadas por √y. En rojo se muestra la función en una escala lineal, mientras que en azul se muestra la función en una escala logarítmica.
Sección 7. Ajuste de cuadrados mínimos a una función arbitraria

Lectura Obligatoria: Bevington (Capítulo 8)
Lectura Sugerida: Numerical Recipes (Secciones 15.5 y 15.6)

7.1 Ajustes no lineales

Los métodos descritos anteriormente para modelar datos están restringidos a funciones que son lineales en los coeficientes a_j:



De esta manera podemos obtener un set de m ecuaciones lineales para las incógnitas el cual puede ser resuelto analíticamente. A continuación consideraremos funciones y(x) no lineales en los coeficientes. Un ejemplo de dicha función es la siguiente:



Esta función es un buen modelo para el decaimiento radioactivo de dos isótopos de corta vida media. La siguiente figura muestra el número de cuentas detectadas en intervalos de 15 segundos durante un lapso de 4 minutos.


Número de cuentas detectadas de una fuente radioactiva en función del tiempo. Se asumen incertezas estadísticas. La curva se obtuvo realizando un ajuste de cuadrados mínimos no lineal.

El modelo considerado tiene 5 coeficientes y contiene tres términos que corresponden a la radiación de fondo (en ausencia de los dos isótopos), la radiación de un isótopo y la radiación del segundo isótopo.

A continuación describimos la generalización del método de máxima probabilidad para el caso no lineal. La probabilidad de observar el set de N observaciones es:



El método consiste en maximizar esta función, o bien minimizar el exponente:



Los valores óptimos de los coeficientes se obtienen cuando las derivadas parciales de χ² se anulan:



Este es un set de m ecuaciones NO lineales en los coeficientes. La manera de resolverlo es considerar a χ² como una función continua en los coeficientes la cual describe una hipersuperficie en un espacio de m dimensiones y buscar en dicho espacio el mínimo de χ². Esta figura muestra la hipersuperficie para el caso de dos parámetros:


Hiper-Superficie de χ² como función de dos parámetros.

7.2 Variaciones de χ² alrededor del mínimo

Para el caso de una muestra grande, la función χ² se comporta como una Gaussiana en cada parámetro alrededor del mínimo.



Vemos que χ² se puede expresar como:



con lo cual se puede ver la variación de χ² con respecto a un parámetro individual cerca del mínimo.



Un alejamiento de una desviación estándar en el parámetro a_j desde el mínimo conlleva un aumento de una unidad en χ².

Este resultado es consistente con aquel obtenido de realizar una expansión de Taylor de χ² alrededor del mínimo:



Como las derivadas parciales alrededor del mínimo se anulan, podemos esperar que χ² sea una función cuadrática cerca del mínimo.

Podemos obtener otra relación útil de la segunda derivada de χ² con respecto al parámetro a_j:



Con esto obtenemos la incerteza en a_j en términos de la curvatura de la función χ² en la región del mínimo:



Basta con tener tres valores de χ² alrededor del mínimo, χ²₁=χ²(a_j1), χ²₂=χ²(a_j2), χ²₃=χ²(a_j3), (en donde a_j2=a_j1+Δa_j y a_j3=a_j2+Δa_j) y realizar un ajuste de una parábola para determinar el mínimo para el coeficiente a_j. La solución está dada por:



Este método se ilustra gráficamente para el caso anterior del decaimiento radioactivo:


Gráfico de χ² como función de un solo parámetro (a₅) alrededor del mínimo. La ubicación del mínimo se calculó ajustando una parábola a los tres puntos.

Además podemos estimar los errores en a_j variando el parámetro hasta que χ² aumente en una unidad desde el valor mínimo.



7.3 Búsqueda del mínimo χ²

Para el caso no lineal existen varios métodos numéricos para encontrar el mínimo χ². Primero debemos partir de valores iniciales razonables y luego calcular χ² en pasos pequeños.

Una complicación es que puede haber más de una solución debido a la presencia de mínimos locales dentro del rango de parámetros considerado:


Hiper-Superficie de χ² como función de dos parámetros.

El método "GRID-SEARCH" consiste en buscar el mínimo separadamente para un parámetro. Cuando el mínimo en esa dirección es hallado se fija su valor y se busca el mínimo en la otra dirección, y así sucesivamente para los demás parámetros. El proceso se repite hasta converger al mínimo.

La siguiente figura muestra un gráfico de contornos para el caso de dos parámetros. Cada contorno muestra niveles de χ² constante. Estos contornos reciben el nombre de contornos de confianza (confidence contours). La línea de zigzag representa el camino de aproximación al mínimo usando el método "GRID-SEARCH".


Gráfico de contorno de χ² como función de dos parámetros (a₁ y a₂). El trazo de zigzag representa el camino de búsqueda del mínimo local unsando el método "GRID-SEARCH".

Este ejemplo muestra que los contornos son altamente elípticos cerca del mínimo. El grado de inclinación de la elipse indica el grado de correlación de los dos parámetros. Si los dos parámetros no estuvieran correlacionados la variación de χ² con cada parámetro debiera ser independiente de la variación del otro parámetro, por lo que el eje de la elipse debiera ser paralelo a cada eje de coordenadas.

Existen otros métodos más eficientes de búsqueda del mínimo tales como el método Lavenberg-Marquardt (ver Numerical Recipes).

Sección 8. Calidad del ajuste

Lectura Obligatoria: Bevington (Capítulo 11)
Lectura Sugerida:

8.1 El test de χ²

La varianza del ajuste s² está dada por:



en donde el factor de peso para cada medición es:



La definición de χ² es:



La relación entre χ² y s² es:



en donde <σ_i²> es el promedio ponderado de las varianzas individuales:



La varianza del ajuste s² caracteriza tanto la dispersión de las mediciones como la calidad del ajuste. La definición de χ², como la razón de s² y de la varianza de la muestra madre σ², es una medida de la calidad del ajuste ("goodness of fit").

Si el modelo de ajuste es una buena descripción de los datos esperamos que la varianza del ajuste s² concuerde con la varianza de la muestra madre σ² de modo que χ² reducido debe ser aproximadamente uno (χ_ν²=1). Si el ajuste es malo, entonces esperamos que χ_ν²>1. Un valor de χ_ν² muy por debajo de 1 indica que hemos sobre-estimado los errores de medición.

La función de distribución de probabilidad para χ_ν² con ν grados de libertad está dada por:



La probabilidad de obtener un valor entre x²=χ² y x²=∞ es:



Esta ecuación da la probabilidad de que un set aleatorio de n mediciones extraídas de la muestra madre arroje un valor de χ² mayor o igual al valor tabulado. Para un modelo de ajuste que se aproxima bien a la función madre el valor de χ_ν² esperado es aproximadamente 0.5 (50%). Para ajustes peores el valor de χ_ν² aumentará y la probabilidad disminuirá. Los valores de la probabilidad integrada se encuentran en la tabla C.4.

8.2 El test F

El test de χ² mide tanto la discrepancia entre los datos propiamente tales como la discrepancia entre el modelo adoptado y la verdadera distribución madre. Un test estadístico que permite separar ambos efectos es el test F.

Si uno determina dos valores χ₁² y χ₂² que siguen la distribución χ², la razón entre χ_ν1² y χ_ν2² se distribuye de acuerdo a la distribución F:



cuya forma funcional es:



La razón siguiente también se rige de acuerdo a la distribución F:



Tal como con el test χ², estamos interesados en la probabilidad integral:



Esta ecuación da la probabilidad de que un set aleatorio de n mediciones extraídas de la muestra madre arroje un valor de F mayor o igual al valor tabulado. Los valores de estas probabilidades se encuentran en la tabla C.5.

Típicamente este test se usa para justificar la inclusión de un término adicional a un ajuste polinomial. Para esto se efectúa un ajuste con m parámetros (N-m grados de libertad) y se calcula χ²(m). Luego agregamos un nuevo parámetro y volvemos a calcular χ²(m+1). La razón Fχ



mide cuanto ha mejorado (disminuido) χ² al agregar un coeficiente. Un valor pequeño para el cambio en χ² sugiere que el término adicional no se justifica. La tabla C.5 da la probablididad de exceder el valor medido de Fχ por azar. Si esta probabilidad es baja (Fχ grande) entonces se justifica agregar el nuevo coeficiente en el ajuste.

8.3 Contornos de confianza

Generalmente asumimos la distribución de Gauss para las incertezas:



Aproximadamente 68% de los eventos se espera que caigan dentro de ±1σ de la media y 95% de las veces dentro de ±2σ. Esto es, podemos afirmar con un 68% de confianza que el valor verdadero esté entre <x>+s y <x>-s.

En general, para una distribución P_x(x;μ) podemos calcular la probabilidad de que una medición caiga entre <x>-a y <x>+b:



y determinar así un grado de confianza para sostener que el valor verdadero de x se encuentre en este rango.

En general la función de probabilidad no tiene porqué ser simétrica por lo que a puede ser diferente de b.

En este gráfico se muestran los contornos correspondientes a 68.3%, 95.4% y 99.73% de confianza para la determinación de los parámetros cosmológicos (Ω_M,Ω_Λ). Los contornos cubren rangos bien diferentes a lo largo de cada eje de coordenadas debido a la alta correlación entre los parámetros. Para este caso en que tenemos dos grados de libertad, la tabla C.4 nos dice que el contorno de 68.3% de confianza corresponde al nivel de chi cuadrado 2.30 unidades mayor que el minimo; para el contorno de 95.4% el nivel con chi cuadrado es 6.17 unidades mayor que el mínimo. Para el contorno de 99.73% el nivel de chi cuadrado es 11.8 unidades mayor que el mínimo. La proyeccion del contorno de confianza del 68.3% en el eje de cada parametro permite especificar el rango del parametro que contiene un 68.3% de probabilidad.

Si pudieramos fijar el valor de uno de los dos parametros (por ejemplo, por medio de un experimento independiente), el valor del otro parametro queda mejor acotado. En este cado tenemos un solo grado de libertad, en cuyo caso los contornos de confianza que contienen probabilidades de 68.3%, 95.4% y 99.73% corresponden a niveles de chi cuadrado 1, 4 y 9 unidades mayores que el minimo, respectivamente.


Gráfico de contornos de χ² como función de dos parámetros (Ω_M,Ω_Λ). Se puede ver un alto grado de correlación entre los dos parámetros. Los contornos corresponden a 68%, 95% y 99% de confianza.

8.4 Test de Kolmogorov-Smirnov (K-S)

Este test se usa para comparar una distribución empírica con una distribución madre. Este test requiere convertir las dos distribuciones a histogramas cumulativos. El test K-S mide la máxima diferencia, D, entre las dos distribuciones cumulativas y evalúa la probabilidad de que esta diferencia ocurra por azar al extraer datos aleatorios de la muestra madre. Este método se puede usar también para comparar dos distribuciones empíricas en que una de las dos distribuciones juega el rol de distribucion madre.
Sección 9. Simulaciones Monte Carlo (EN CONSTRUCCION)

Lectura Obligatoria: Bevington (Capítulo 5)
Lectura Sugerida:

Generalmente, los experimentos son el resultado de integraciones de multiples distribuciones de probabilidad que no tienen una evaluacion analitica. Una posibilidad es integrar las integrales numericamente evaluando la funcion sobre una grilla regular de N elementos. Cuando hay m variables involucradas, la dificultad numerica aumenta como N^m. Las simulaciones Monte Carlo se usan para evaluar integrales usando un muestreo aleatorio de las variables en vez de una grilla regular.

Por ejemplo, si estamos interesados en averiguar el area de un circulo, podemos generar 100 eventos aleatorios entre x[-1,+1] y[-1,+1] y calcular la fraccion de eventos que cayeron dentro del circulo. Esta simulacion la podemos repetir 100 veces y construir el siguiente histograma, con lo cual podemos estimar trivialmente la incerteza en el area del circulo.


Histograma del área de un círculo obtenida a partir de 100 simulaciones independientes Monte Carlo, cada una con 100 pares de números aleatorios. En azul se muestra una distribución Gaussiana con μ=3.127 y σ=0.156 obtenida a partir de las 100 simulaciones.

Los computadores estan dotados de rutina que generan numeros aleatorios ("uniform deviates") en el rango r=0-1. Si definimos p(r) como una distribucion de densidad de probabilidad uniforme entre 0 y 1:



La distribucion esta normalizada:



Si quisieramos generar una distribucion uniforme entre x=-1 y x=1,



necesitaremos calcular x a partir de r:



Este es un ejemplo simple de transformacion lineal de variables.

En general estaremos interesados en generar una numeros aleatorios para distribuciones no uniformes (e.g. una distribucion normal). Esto se puede lograr haciendo transformacion de variables a partir de r, teniendo en cuenta que la probabilidad integrada en cualquiera de ambas variables debe ser igual.



Por tanto, podemos escribir



Usando el hecho que P(r)=1 entre r=0 y 1:



Dado r (generado aleatoriamente a partir de la distribucion uniforme) y la nueva distribucion de probabilidad que nos interesa, P(x), podemos calcular numericamente esta integral y determinar el valor de x que satisface esta ecuacion (a veces la integral admite una solucion analitica).



El valor de x es un numero aleatorio extraido de una distribucion no uniforme P(x) por el "metodo de transformacion de variables".

Este metodo se puede aplicar, por ejemplo, para la distribucion Gausiana y generar desviaciones gausianas.



Tambien se puede aplicar a la distribucion de Poisson



reemplazando la integral por una suma. Este histograma fue generado por la tecnica de Monte Carlo para una distribucion de Poisson con una media de μ=8.4.


Histograma de 200 variables aleatorias generadas de una distribución Poissoniana con μ=8.4.

A menudo se requiere extraer numeros de una distribucion exponencial:



El metodo de transformacion de variables arroja una expresion analitica para la variable t a partir de numeros aleatorios r extraidos de la distribucion uniforme:



mhamuy@das.uchile.cl