Banco Central de Venezuela

Colección Economía y Finanzas Serie Documentos de Trabajo

SELECCIÓN DE VARIABLES POR MUESTREO PARA REGRESIÓN MÚLTIPLE BASADO EN EL RISK INFLATION CRITERION MIGUEL DORTA

[Nº 77] Mayo 2006

Banco Central de Venezuela

Colección Economía y Finanzas Serie Documentos de Trabajo

SELECCIÓN DE VARIABLES POR MUESTREO PARA REGRESIÓN MÚLTIPLE BASADO EN EL RISK INFLATION CRITERION MIGUEL DORTA [email protected]

[Nº 77] Mayo 2006

© Banco Central de Venezuela, Caracas, 2008 Gerencia de Investigaciones Económicas Producción editorial Gerencia de Comunicaciones Institucionales Departamento de Publicaciones Torre Financiera, piso 14, ala sur. Avenida Urdaneta, esquina de Las Carmelitas Caracas 1010 Teléfonos: 801.8075 / 8063 Fax: 536.9357 [email protected] www.bcv.org.ve

Las opiniones y análisis que aparecen en la Serie Documentos de Trabajo son responsabilidad de los autores y no necesariamente coinciden con las del Banco Central de Venezuela.

Se permite la reproducción parcial o total siempre que se mencione la fuente y no se modifique la información.

ÍNDICE

1.

INTRODUCCIÓN................................................................................ 9

2.

EL RISK INFLATION CRITERION PARA REGRESIÓN MÚLTIPLE .......................................................11

3.

DESCRIPCIÓN DEL ALGORITMO SVMRIC ............................... 13 3.1 El algoritmo .................................................................................. 13 3.2 Definiciones .................................................................................. 13 3.3 Los pasos del algoritmo ................................................................ 14

4.

EXPERIMENTO MONTE CARLO PARA LA ESTIMACIÓN DE UNA FUNCIÓN DE LOS PARÁMETROS INICIALES......................................................................................... 17

5.

EJEMPLOS ........................................................................................ 22 5.1 Ejemplo: No colinealidad, R2 ajustado alrededor de 0.86 y T=100 .............................................................................................. 23 5.2 Ejemplo: Combinaciones de varios niveles de colinealidad, R2 , tamaños de muestra y número de iteraciones .............................. 24

6.

CONCLUSIONES.............................................................................. 27 APÉNDICE A..................................................................................... 28 APÉNDICE B..................................................................................... 30 REFERENCIAS BIBLIOGRÁFICAS ............................................... 34

Resumen En este trabajo se desarrolló un algoritmo que permite seleccionar y estimar una ecuación dinámica de regresión lineal por OLS a partir de un número grande de posibles términos explicativos. El algoritmo va seleccionando por muestreo Bernoulli con reemplazamiento los términos explicativos con probabilidades que se actualizan progresivamente en función de si cumplen o no el Risk Inflation Criterion (Foster y George, 1994), evitando así la búsqueda exhaustiva entre todos los subconjuntos posibles de términos. Se desarrollan experimentos Monte Carlo tanto para obtener una fórmula empírica para la inicialización de los parámetros del algoritmo adecuados a diferentes casos y también para producir ejemplos, que permitan tener una idea del desempeño del algoritmo bajo diferentes condiciones.

Palabras clave: Selección Estadística de Modelo, Selección de Variables, Modelo de Regresión Lineal, Criterios Estadísticos de Información, Simulaciones Monte Carlo, Muestreo. Clasificación JEL: C15, C22, C51

Abstract This document reports the development of an algorithm that allows the selection and estimation by OLS of a linear dynamic equation departing from a large number of possible explicative terms. The algorithm sequentially selects the explicative terms by Bernoulli sampling with dynamically adjusted probabilities depending on whether the Risk Inflation Criterion (Foster y George, 1994) is satisfied, avoiding in this way the exhaustive search among all possible subsets of explicative terms. Monte Carlo experiments are design and run in order to obtain an empirical formula to give initial values to the algorithm’s parameters adequate to a range of cases and also for producing examples that permit to have an idea of the algorithm’s performance under different conditions. Key words: Statistical Model Selection, Variable Selection, Linear Regression Model, Statistical Information Criteria, Monte Carlo Simulations, Sampling. JEL classification: C15, C22, C51

1. INTRODUCCIÓN

El problema de selección estadística de variables para un modelo de regresión lineal es de tal relevancia que ha propiciado la creación de una impresionante variedad de métodos. Miller (1990), recopila y trata excelentemente los principales métodos previos a ese año. Uno de los principales aspectos de la selección de variables es su enorme tamaño. Aún con moderadas cantidades de variables, el análisis del espacio de todas las regresiones posibles demanda demasiado tiempo de procesamiento dado el presente estado del arte en computación1. En este sentido, muchas estrategias se basan en la reducción de dicho espacio como los métodos de Efromyson (1966) y Furnival y Wilson (1974). Posteriormente, sin dejar de prestar atención a la reducción del espacio de modelos, se comienzan a profundizar y mejorar en los criterios estadísticos de selección. Los más familiares quizás son los propuestos por Mallows (1973), Akaike (1973) y Schwarz (1978), los cuales se basan en la familia de criterios de la suma de cuadrados penalizada. Más recientemente se han propuesto nuevos ajustes basados en esta familia de criterios, entre los cuales figura el Risk Inflation Criterion (RIC) propuesto por Foster y George (1994) y que es fundamental en el algoritmo que se propone en el presente trabajo. También vinculados a la misma familia existen otros métodos y criterios de selección como los sugeridos por Tibshirani y Knight (1999), Ye (1998), Hurvitz y Tsai (1989, 1998), Rao y Wu (1989), Shao (1997), Wei (1992), Zheng y Loh (1997), Benjamini y Hochberg (1995), Clyde y George (1999, 2000) y George y Foster (2000). Si nos extendiéramos a otras estrategias de selección estadística de variables para modelos de regresión lineal, las referencias serían mucho más extensas, pero no lo hacemos debido a que ellas están menos relacionadas al algoritmo que se desarrolló en este trabajo.

1

Por ejemplo, para un problema relativamente pequeño como de una población de 30 posibles términos explicativos habría que correr 230 = 1.073.741.824 regresiones.

9

Para el caso que se desarrolla en esta investigación comenzamos por suponer que los datos a ser analizados son series de tiempo con las cuales se desea estimar una ecuación dinámica de regresión lineal por OLS de una variable dependiente para la cual se dispone de una “población” de potenciales términos explicativos (variables explicativas y sus rezagos). Supongamos, además, que el número de términos explicativos es considerablemente grande (incluso mayor que el tamaño de la muestra) y que puede contener muchos términos irrelevantes. Una forma de proceder es efectuar una búsqueda exhaustiva entre todos los subconjuntos posibles y seleccionar aquel que cumpla con algún criterio estadístico de información. Pero en un caso donde la “población” de términos es muy grande, el tiempo de procesamiento puede ser prohibitivamente costoso. Con este problema en mente, una manera intuitiva de proceder es diseñar un algoritmo que trabaje de manera similar a lo que haría el cerebro humano. Es decir, si el investigador no tiene ninguna información a priori sobre cuáles son los términos relevantes y cuáles no, probablemente comenzaría a realizar una suerte de muestreo donde las probabilidades de selección iniciales de las variables van a ser ajustadas de acuerdo con la significación estadística que va observando en la secuencia de muestreo y ajustes por OLS. Entonces, el investigador tenderá a seleccionar más frecuentemente los términos que en repetidas ocasiones mantienen su significación, alegando robustez estadística y probará con menor frecuencia los términos que reportan en reiteradas oportunidades menor significación estadística. En ese proceso, cada vez se le irá dificultando el poder superar el anterior mejor modelo y eventualmente parará su búsqueda, considerando que ese modelo es apropiado para explicar el fenómeno en estudio2. Por otra parte, al ser la “población” de términos explicativos muy grande, el investigador corre el riesgo de seleccionar en exceso términos irrelevantes, aún con niveles de significación típicamente usados como 5% e incluso 1%. Por ejemplo, si la “población” de términos explicativos es de tamaño 200 y todas a priori se conocen como irre2

Un algoritmo que presenta algunas similitudes al propuesto en este trabajo puede ser encontrado en George y McCulloch (1993).

10

levantes para explicar la variable dependiente, el uso de un nivel de significación de 5% tendería a seleccionar modelos con un promedio de 10 variables explicativas, cuando en realidad no debería seleccionar variable alguna. En este orden de ideas, el procedimiento de muestreo mencionado parece ser idóneo si es combinado con el criterio RIC. Foster y George (1994), desde un enfoque minimax de la teoría de la decisión, demuestran que el criterio RIC permite la más pequeña posible inflación máxima en riesgo predictivo debida a selección, cuando la dimensión de la población de términos se hace muy grande. El resto del trabajo se organiza de la manera siguiente: En la segunda sección se explica en más detalles el criterio RIC. En la tercera sección se describe el algoritmo de selección de variables por muestreo para regresión múltiple basado en ese criterio (SVMRIC). En la cuarta sección, se describe el diseño y los resultados de un experimento Monte Carlo (MC) sobre el algoritmo SVMRIC que tiene por finalidad obtener una ecuación general para la estimación de parámetros iniciales en función de características de los datos y del criterio de parada del algoritmo. En la quinta sección se desarrollan ejemplos simulados que permiten obtener una idea del desempeño del algoritmo bajo diferentes condiciones. En la última sección se ofrecen algunas conclusiones generales y sugerencias tanto para el uso del algoritmo como para futuras investigaciones. 2.

EL RISK INFLATION CRITERION PARA REGRESIÓN MÚLTIPLE

En virtud de que el criterio RIC es fundamental para este trabajo, en esta sección se trata de formar una intuición sobre el mismo de manera similar a como lo explica George (2000). Para el modelo lineal, muchos de los criterios más populares de selección son casos especiales del criterio de la suma de cuadrados penalizada, el cual constituye un marco unificado de comparación.

11

Asumiendo que σ es conocido para evitar complicaciones, este criterio general selecciona el subconjunto de variables que minimiza (RSSγ/ σ2 + Fqγ), donde F es una dimensión de penalización preestablecida. Esta expresión penaliza RSSγ/σ2 por F veces qγ , la dimensión del γ-esimo modelo. El criterio AIC corresponde a F = 2, y el criterio BIC se obtiene fijando F = log T, donde T es el tamaño de muestra. A menos que T sea muy pequeño, AIC impone una penalización más chica y por lo tanto seleccionará modelos más grandes que BIC. Si suponemos ahora que hay p predictores ortogonales, entonces la suma de cuadrados penalizada selecciona aquellos predictores con t-estadísticos tal que t2>F. Si hay muchos predictores y estos no están relacionados con la variable dependiente (todos los p coeficientes de la regresión son 0) entonces AIC es claramente muy liberal y tenderá a incluir muchos predictores irrelevantes. Bajo este supuesto, cuando p es grande el valor esperado de t2 es aproximadamente 2logp. En este sentido, una selección natural más conservativa es F = 2logp. Este es el criterio RIC propuesto por Foster y George (1994) y también es el umbral universal para ondículas (wavelets) propuesto por Donoho y Johnstone (1994). Estos dos artículos motivaron, desde un enfoque minimax de la teoría de la decisión, a F = 2logp como el que da la más pequeña posible inflación máxima en riesgo predictivo debida a selección (cuando p→∞). Note que BIC coincide con RIC si 2logp=logT y por lo tanto también será muy liberal si p2>T. Para p muy grande, lo anteriormente expuesto implica que los niveles de significación usados para descartar términos irrelevantes deberían en muchos casos ser más pequeños que los usados comúnmente en la estadística aplicada. Por ejemplo: si el período muestral es suficientemente grande como para usar la aproximación normal de la distribución t, ya para p=7 el criterio RIC implica un nivel de significación de 4,85% que es menor al típico 5%, para p=28 implica un nivel de significación de 0,98% menor al típico 1%. En un caso donde p=100 el nivel de significación implicado por el RIC sería de apenas 0,24%.

12

3.

DESCRIPCIÓN DEL ALGORITMO SVMRIC

El algoritmo que aquí proponemos lo que pretende hacer es combinar las bondades del criterio RIC con un procedimiento que permita reducir sensiblemente el espacio de modelos a ser comparados con dicho criterio. Omitiendo algunos detalles, dicho en palabras lo que el procedimiento propuesto hace en esencia es realizar una secuencia de selección de subconjuntos de términos explicativos por muestreo para ser probados en regresión, tal que la probabilidad de selección inicial de cada término explicativo en la población va a ser ajustada de acuerdo con el número de veces en que cada término cumple con el criterio RIC. Entonces, se tenderá a seleccionar cada vez con mayor frecuencia a los términos que más cumplan con ese criterio en las pruebas sucesivas acumuladas, y con menor frecuencia a los términos que reportan en reiteradas oportunidades menor significación estadística relativa al RIC. En ese proceso, las frecuencias relativas de cumplimiento del criterio RIC se estabilizan y con ello se tiende a converger a un modelo que tiene un máximo de predictores que cumplen con el RIC. 3.1 El algoritmo

Se desea estimar una ecuación de regresión por OLS para una variable dependiente Yt sobre la cual existe una “población” de potenciales términos explicativos {Xi,t-j} donde i=1,…,q es el subíndice referido a las posibles variables explicativas, j=0,1,…,r es el subíndice referido a los posibles rezagos y t=1,…,T es el período muestral3. 3.2 Definiciones

Sean los elementos de una matriz (P) que en la iteración h contiene las probabilidades de seleccionar el término Xi,t-j. 3

La matriz X puede o no contener rezagos de la variable dependiente y términos contemporáneos de otras variables dependiendo de si estas son o no exógenas.

13

Sean los elementos de una matriz (F) que en la iteración h contiene el número de veces (o frecuencia absoluta) en que el término Xi,t-j ha sido seleccionado en la muestra para ser probado en la regresión. Sean los elementos de una matriz (G) que en la iteración h contiene el número de veces con una pequeña corrección (frecuencia absoluta ajustada) en que un término Xi,t-j, seleccionado para ser probado en la regresión, cumpla con el criterio RIC. Sea p la probabilidad inicial de seleccionar un término para ser incluido en una regresión. Sea pl un parámetro que, como veremos más adelante, permite regular las probabilidades mínimas de que un término sea seleccionado en la regresión, cuando el algoritmo se acerca hacia la convergencia de P. 3.3 Los pasos del algoritmo

=p,

Se inicializan Note que

/

=q.r,

=p.q.r para todo i, j

=p para todo i, j

Paso h

Se genera la matriz cuyos elementos contienen realizaciones de variables aleatorias independientes Bernoulli de parámetros . Entonces los elementos de esta matriz contienen unos o ceros lo que define si el término explicativo fue o no seleccionado para regresión, respectivamente.

14

Luego se ajusta por OLS4 la ecuación (1) que no es otra cosa que la regresión con solo la muestra de términos seleccionados. Paso h+1

Si entonces y si el coeficiente del término 2 correspondiente cumple con que t >2*log(q.r), es decir, si se cumple el criterio RIC, calcular . En caso de no cumplir el criterio RIC, calcular , donde es el p-valor que mide la significación estadística del término correspondiente basa. do en su estadístico t5. Luego se calcula Si alguna de las variables seleccionadas no cumple con el criterio RIC, se vuelve al paso h sustituyendo h por h+1. En caso contrario, se compara el coeficiente R2 ajustado con el de la anterior mejor regresión y se registra como mejor regresión a la que tenga el mayor coeficiente6. Si no se encuentra una regresión mejor en este paso, entonces se calcula kh+1=kh+1 y se vuelve al paso h sustituyendo h por h+1. En el caso contrario se reinicializa kh+1=0 y se vuelve al paso h sustituyendo 4 5

6

Parece recomendable usar errores estándares robustos a heterocedasticidad y autocorrelación según Newey-West. Sin embargo, hay que tener en cuenta que los tiempos de procesamiento se pueden incrementar considerablemente. La intuición es que, hacia la convergencia de la matriz P, la probabilidad de seleccionar cualquier variable no pueda caer por debajo de un mínimo que depende de pl y de la secuencia de significaciones estadísticas de los términos de las regresiones que no cumplan el criterio RIC. Note que si no se efectuara este cálculo, los elementos de G serían frecuencias absolutas y de allí que las hemos denominado “frecuencias absolutas ajustadas”. El óptimo de p y pl dependería de la estructura de los datos y del criterio de parada del algoritmo. Esto se investiga con experimentos MC los cuales se presentan más adelante. Note que se está maximizando el R2 ajustado dentro de la clase de regresiones cuyos términos cumplen individualmente con el criterio RIC. Para q.r>T el resultado combinado será más conservador que el uso aislado del criterio de Schwarz.

15

h por h+1. El algoritmo se detiene luego de que kh exceda un número predeterminado n que entonces pudiera denominarse como el “número de iteraciones sin mejora”. Inmediatamente surge la pregunta: ¿Cómo fijar los parámetros iniciales n, p y pl? Como veremos más adelante la respuesta empírica a esa pregunta se intenta resolver con un experimento MC. Sin embargo, previamente vale la pena desarrollar algunas intuiciones. En principio, n debería ser lo más grande posible y en este sentido, esto debería obedecer a razones prácticas de costo en el tiempo de procesamiento que la aplicación particular pueda permitir por razones técnicas y prácticas. En las pruebas MC, que se describen en la próxima sección, el algoritmo mostró una velocidad de convergencia relativamente rápida de la matriz P. Sin embargo, una vez que la matriz converge no necesariamente ya se ha encontrado una regresión “óptima” o muy cercana a ella. Es por ello que el criterio de parada del algoritmo se considera como un número predeterminado n de iteraciones sin encontrar una regresión superior a la anterior mejor regresión. La probabilidad inicial de selección de un término (p) determina la dimensión de los primeros modelos de regresión probados por el algoritmo. Así, valores grandes de p implican que las primeras regresiones tendrían muchos términos, en tanto que valores pequeños de p conducen a que las primeras regresiones sean de pocos términos. Por otra parte, el parámetro pl, modula las probabilidades mínimas (hacia la convergencia de la matriz P) de que un término sea seleccionado para ser incluido en una regresión. La magnitud “óptima” o adecuada de estos parámetros debería depender de características de los datos para lo cual se emplean los autovalores de la matriz de correlación de la población de términos potencialmente explicativos. Finalmente, vale la pena comentar que dado un valor de n, posiblemente los parámetros iniciales “óptimos” p y pl presenten mutua dependencia.

16

Esto es importante a la hora de tratar de estimar alguna función que permita aproximar dichos “óptimos” y para lo cual en este trabajo se utiliza un diseño de experimento MC. 4. EXPERIMENTO MONTE CARLO PARA LA ESTIMACIÓN DE UNA FUNCIÓN DE LOS PARÁMETROS INICIALES

Este experimento tiene como principal objetivo el poder dar estadísticamente con una función lo más general posible para la estimación de los parámetros iniciales “óptimos” o adecuados del algoritmo SVMRIC que dependa de las características de los datos. Por razones prácticas es necesario que esta estimación de parámetros iniciales sea muy rápida en relación al tiempo medio que se espera pueda consumir el algoritmo SVMRIC en su ejecución. Así el conjunto de variables que potencialmente puedan ayudar a obtener esta función también debe ser de muy rápida obtención. El experimento consiste en tratar de obtener una matriz (XXP) suficientemente grande que contenga parámetros p y pl “óptimos” acompañados de un conjunto de variables que resuman características de la data y que potencialmente pueden explicar a dichos parámetros. Es deseable que el algoritmo SVMRIC disponga inicialmente de una función general de estimación rápida de los parámetros p y pl que puedan ser desde el punto de vista práctico suficientemente “buenos” para aproximarse con mayor probabilidad al modelo verdadero, para una gama muy variada de procesos generadores de datos (PGD). En este sentido, se requiere que cada observación generada para la matriz XXP provenga de PGDs de características muy variantes. Para esto, se diseñó una subrutina de generación de datos (GENDATA7) que produce, de una corrida a otra, datos y modelos con diferentes estructuras. Así, para producir una observación de la matriz XXP, se comienza por correr la rutina GENDATA después de haber generado aleatoriamente 7

La descripción detallada de esta subrutina se encuentra en el Apéndice A.

17

una serie de parámetros que la condicionan como lo son: el tamaño de muestra (tm) que es extraído de la distribución Ud(100,300)8, el tamaño “poblacional” de variables explicativas (nzmx) generado de la distribución Ud(6,9), el tamaño “poblacional” de rezagos (nrzmx) extraído de la distribución Ud(6,18), la desviación estándar de la regresión verdadera (sx1) que se genera de la distribución Uc(0,2), el número de variables de la regresión verdadera (nz) extraído de la distribución Ud(1,nzmx) y el número de rezagos de la regresión verdadera (nrz) extraído de distribución Ud(1, entero(0.5nrzmx)+1). Fijos estos datos y parámetros, se genera un n de la distribución Ud(500,5500) y se realizan un número suficientemente grande de replicaciones del algoritmo SVMRIC9 con valores de p y pl generados aleatoriamente10. De este conjunto de replicaciones se identifica cuál es la que produce el mejor modelo de regresión11 llenando así una observación (fila) de la matriz XXP. Una vez terminadas las simulaciones MC, queda todo dispuesto para proceder a estimar las funciones de inicialización de p y de pl. Las columnas de la matriz XXP (1200 filas u observaciones) están conformadas por variables que pueden clasificarse en dos tipos: 8 9

10

11

Ud y Uc denotan a las distribuciones uniformes discreta y continua, respectivamente. En este trabajo el número de replicaciones varió entre un mínimo de 30 y un máximo de 100, dependiendo de si en la secuencia no se superaba el último mejor modelo luego de 30 replicaciones. Aunque, en un principio, estos números no parecieran lo bastante grandes, el modo en que durante el desarrollo del experimento se iba ajustando la distribución de donde se extraían los valores de p y pl, los hacen suficientes desde el punto de vista práctico. Adicionalmente, ya con esta configuración se alcanzaba el tiempo de procesamiento previsto para este estudio. Las primeras 20 observaciones de p y pl en la matriz XXP se obtienen a partir de distribuciones uniformes U(0,1) usadas en las replicaciones del algoritmo SVMRIC. No obstante, note que los “óptimos” encontrados en cada conjunto de replicaciones provienen de distribuciones desconocidas que dependen en gran medida de la subrutina GENDATA y de la previa generación de parámetros que la afectan. Así, la generación de p y pl a partir de la observación 20 se realiza con distribuciones Beta cuyos parámetros se estiman por el método de los momentos siendo actualizados en cada nueva observación generada para XXP. Como es de esperar, son necesarios cierto número de generación de observaciones de la matriz XXP para que los parámetros de la distribución Beta se estabilicen adecuadamente. Así, en este trabajo se desecharon alrededor de las primeras 80 observaciones para asegurar una estabilización adecuada. Dentro de la clase de modelos donde todos sus términos explicativos cumplen con el criterio RIC, se toma el que tenga mayor R2 ajustado y que se encuentre con el menor número de iteraciones.

18

Variables posiblemente explicativas de p y/o pl

tm: Tamaño de muestra de la data generada por la subrutina GENDATA n: Número máximo de iteraciones sin mejora del algoritmo SVMRIC p: Probabilidad inicial de que una variable sea seleccionada. pl: Parámetro que modula la probabilidad mínima de que una variable sea seleccionada hacia la convergencia del algoritmo SVMRIC mean: Promedio de los autovalores (en proporción) de la matriz de correlación de los datos generados por la subrutina GENDATA median: Mediana de los autovalores (en proporción) de la matriz de correlación de los datos generados por la subrutina GENDATA max: Máximo de los autovalores (en proporción) de la matriz de correlación de los datos generados por la subrutina GENDATA min: Mínimo de los autovalores (en proporción) de la matriz de correlación de los datos generados por la subrutina GENDATA stdv: Desviación estándar de los autovalores (en proporción) de la matriz de correlación de los datos generados por la subrutina GENDATA nzmx: Tamaño poblacional de variables explicativas nrzmx: Tamaño poblacional de rezagos

19

Otras Variables que se guardaron en la matriz XXP

jk: Número de replicaciones del algoritmo SVMRIC sx1: Desviación estándar de la ecuación de la regresión verdadera nz: Número de variables de la ecuación de la regresión verdadera nzr: Número de rezagos de la ecuación de la regresión verdadera ncoefz: Número de coeficientes de la mejor especificación por algoritmo SVMRIC controlz: Número de iteraciones del algoritmo SVMRIC al encontrar la última mejor regresión rbar2z: R2 ajustado de la mejor ecuación por selección algoritmo SVMRIC dcoefz: Proporción de coeficientes seleccionados respecto de verdaderos Para especificar y estimar las regresiones de p y pl, es conveniente considerar que los valores que pueden tomar estos parámetros deben estar restringidos, por lo menos, al intervalo (0,1). Una especificación apropiada al caso es la siguiente: (2) Esta expresión se puede linealizar mediante la transformación logística lo que conduce a la especificación: (3)

20

Esta misma especificación se usa para p y pl. En lo sucesivo utilizamos la siguiente notación: (4) En este caso, los términos zi son los logaritmos de la lista de variables potencialmente explicativas nombradas anteriormente con la particularidad de que en la ecuación de lpl se incluye adicionalmente a lp como variable explicativa. Esto permite estimar los parámetros p y pl en forma recursiva en las aplicaciones del algoritmo SVMRIC y así resolver de manera práctica el problema de simultaneidad entre p y pl. El producto de nzmax por nzrmax es el número poblacional de términos que es el inverso del promedio de los autovalores de la correspondiente matriz de correlación; y en consecuencia, es preciso excluir indistintamente una de estas dos variables. La estrategia para llegar a las especificaciones y estimaciones finales de las ecuaciones de lp y lpl consiste en comenzar con una especificación general que incluya a todas las variables para luego ir reduciéndola, excluyendo un término a la vez (el menos relevante) hasta que en la ecuación sólo queden términos estadísticamente significativos (también según el criterio RIC). Las estimaciones finales de las ecuaciones para lp y lpl se presentan en el cuadro 1. En la estimación de la ecuación de lp, todas las variables explicativas son significativas a menos del 1% excepto log(n) que lo es al 3%. Sin embargo, el R2 ajustado parece bastante bajo. Por su parte, la estimación de la ecuación de lpl, contiene menos variables explicativas aunque también significativas a menos del 1%. En este caso, el coeficiente R2 ajustado es todavía más bajo que el anterior. Es posible que la naturaleza del problema imponga límites a estos coeficientes de determinación; sin embargo, no se puede descartar que estas estimaciones puedan mejorarse mediante un mejor diseño del experimento MC y/o mediante un análisis más profundo en cuanto a la

21

posible complejidad de la no linealidad de las ecuaciones. Hasta tanto no se logren nuevas y mejores estimaciones, se recomiendan las aquí obtenidas para incorporarlas a programas con fines aplicados. CUADRO 1

5.

EJEMPLOS

En esta sección se intenta demostrar el desempeño del algoritmo SVMRIC mediante ejemplos aplicados sobre datos simulados. Los datos simulados tienen la ventaja de que al conocer a priori el modelo verdadero y la estructura de correlación de las potenciales variables explicativas, se hace posible observar mejor ciertas propiedades del algoritmo. En este sentido, trataremos de estudiar como dichas estructuras de correlación afectan su desempeño así como observar empíricamente aspectos de su comportamiento asintótico. El PGD genérico utilizado para los ejemplos es el siguiente:

22

S0 = Desviación típica muestral de z0,t

El diseño particular de este PGD obedece a la necesidad de variar la estructura de correlación y el coeficiente de determinación de la regresión verdadera. Por ejemplo, para a=0, cov(zi,zj)=0 para todo i≠j ; para a>1,2 se producen estructuras que pudieran considerarse de extrema colinealidad entre los zi . Por otra parte, el R2 de la especificación verdadera dependerá del valor de k. Así, valores de k cercanos a cero (0) producen valores de R2 cercanos a 1, para k=1 el valor de R2 tiende a 0,5; y para k >2 el coeficiente R2 tenderá a cero (0). 5.1 Ejemplo: No colinealidad, R2 ajustado alrededor de 0,86 y T=100

En este ejemplo se usan los siguientes parámetros para el PGD: a=0, k=0,4, T=100, b1=b2=1, b3=-1, b4=0,5, b5=-0,5, b6=0,7, b7=-0,7, b8=0,9 y b9=-0,9. Con fines informativos, en el cuadro 2 se muestran los resultados de los ajustes del modelo verdadero para un par de conjuntos de datos generados con estos parámetros. Luego se generan 100 conjuntos de datos con este PGD y para cada uno de ellos se corre el algoritmo SVMRIC utilizando como población de términos explicativos hasta 8 rezagos de las 9 variables (q.r=72). Habiendo fijado en 500 y 2.000 el número de iteraciones sin mejora (n), el algoritmo acertó la especificación verdadera en 68 y 80 oportunidades, respectivamente. Los promedios de los R2 ajustados para estos aciertos fueron de 0,864 y 0,865, respectivamente. Por su parte, los R2 ajustados promedios de estos coeficientes para todos los modelos seleccionados por el algoritmo fueron, respectivamente 0,861 y 0,863, respectivamente. El algoritmo falló en encontrar solución en solo un caso que correspondió, como es de esperar, al de menos iteraciones.

23

5.2 Ejemplo: Combinaciones de varios niveles de colinealidad, R2, tamaños de muestra y número de iteraciones

En este caso, el modelo verdadero es exactamente igual al del ejemplo anterior (b1=b2=1, b3=-1, b4=0,5, b5=-0,5, b6=0,7, b7=-0,7, b8=0,9 y b9=-0,9). La dimensión de la población de términos explicativos será q=9, r=20 (q.r=180). Se efectúan 100 repeticiones del algoritmo para cada una de las combinaciones provenientes de cruzar: 1) valores de a: {0,0,25,…,1,25}, 2) valores del R2 ajustado configurados por los valores de k: {0,25,0,50,…,1,75}, 3) tamaños de muestra T: {100 y 1.000} y 4) número de iteraciones n: {1.000 y 5.000}. Estos cruces (en total 168) recogen una gama de niveles desde no existencia de colinealidad hasta niveles extremos de la misma, coeficientes R2 ajustados desde muy bajos hasta muy altos y tamaños de muestra con número de iteraciones que tienen el propósito de observar empíricamente el desempeño asintótico del algoritmo. Para cada uno de estos cruces se tabulan las estadísticas siguientes: la frecuencia de aciertos con el modelo verdadero y el correspondiente R2 ajustado promedio; la frecuencia en que el algoritmo encuentra solución y también el R2 ajustado promedio; el R2 promedio del modelo verdadero en todas las replicaciones; el promedio de los máximos de los autovalores; y el promedio de las desviaciones estándar de los autovalores (ver Apéndice B). CUADRO 2

24

Antes de comentar los resultados, es importante destacar la alta complejidad del problema inducida por el alto número (180) de potenciales términos explicativos. El lector puede observar los resultados en detalle en los cuatro (4) cuadros del Apéndice B. Sin embargo, esa información se resume en el cuadro 3 agrupándola, para cada tamaño de muestra, en sextiles según el número de aciertos (dispuestos en orden descendente) y presentando los promedios de las variables. Además, se agregó la diferencia absoluta de los R2 ajustados entre las soluciones encontradas por el algoritmo y los modelos verdaderos. CUADRO 3

Varios aspectos resaltantes se pueden observar en el cuadro 3. En primer lugar destaca el impacto de subir el tamaño de muestra de 100 a 1.000. Se nota como este incremento influye fuertemente en el número promedio de aciertos con el modelo verdadero, sobre todo al comparar los sextiles 1, 2, 3 y 4. Esto parece indicar un buen desempeño asintótico del algoritmo. Además, es importante señalar que las diferencias absolutas de los R2 ajustados promedios entre los casos donde se encuentran soluciones que cumplan los criterios y los modelos verdaderos disminuyen contundentemente cuando se aumenta la muestra a 1.000. En este sentido hay que tomar en cuenta que en los casos donde no se acierta el modelo verdadero, bien se puede haber

25

obtenido un modelo cercano al mismo lo cual puede aproximarse mediante estas diferencias de R2 ajustados. En segundo lugar, al concentrar el análisis para el tamaño de muestra de 100, se observa en el primer sextil un número promedio de aciertos bastante alto (58), pero luego decae dramáticamente a cuatro (4) para el segundo sextil y a cero (0) para el resto de los sextiles. En este sentido, observando los valores promedios de a, de los R2 ajustados del modelo verdadero y de las estadísticas de los autovalores, parece obvio que la mayor colinealidad y/o menores coeficientes de determinación son los causantes de este bajo desempeño a partir del segundo sextil. No obstante, el último sextil reporta baja colinealidad promedio pero en cambio es al que le corresponde el menor número de iteraciones sin mejora promedio y también el penúltimo R2 ajustado promedio de los modelos verdaderos en orden descendente. En tercera instancia, si se concentra el análisis para el tamaño de muestra de 1.000, se aprecian también estos mismos patrones de comportamiento, sólo que en este caso las desmejoras en el desempeño del algoritmo son muy leves en los cuatro primeros sextiles. El bajo desempeño ocurrido en los dos últimos sextiles no sorprende, ya que es donde existen niveles de colinealidad extremos y/o bajo número de iteraciones. Finalmente, parece recomendable el asegurar que el algoritmo se corra con un número suficientemente grande de iteraciones sin mejora para poder obtener resultados satisfactorios. En este sentido, se propone que el investigador utilice como una primera convención informal a n ≥ 1500+10(q.r)12.

12

Estas cifras de iteraciones son ínfimas si se comparan con el número de todos los subconjuntos posibles que para este ejemplo es de 2180 = 1,5(10)54. Para un problema mucho más simple, como el de una población de términos explicativos de tamaño 30, el número de todos los subconjuntos posibles es 230 = 1.073.741.824 que todavía parece muy costoso en términos del tiempo de procesamiento para un computador personal en el presente estado del arte.

26

6.

CONCLUSIONES

En este trabajo se desarrolló un algoritmo que le permite a un investigador seleccionar y estimar una ecuación dinámica de regresión lineal por OLS a partir de un número grande de posibles términos explicativos. El algoritmo va seleccionando por muestreo Bernoulli con reemplazamiento los términos con probabilidades que se actualizan progresivamente en función de si se cumple o no el Risk Inflation Criterion evitando así la búsqueda exhaustiva entre todos los subconjuntos posibles. Después de efectuar un largo y amplio experimento MC para analizar de qué depende la magnitud adecuada de los parámetros iniciales del algoritmo según el tipo de problema, se encontró que los factores relevantes fueron: el tamaño de la “población” de términos explicativos, el máximo de los autovalores (en proporción) de la matriz de correlación de dichos términos, la desviación estándar de esos autovalores y el número de iteraciones sin mejora que se le imponga al algoritmo como criterio de parada. Sin embargo, es importante mencionar que esta relación, aunque robustamente significativa, produce unos R2 ajustados aparentemente bajos. En este sentido, parece haber margen para mejorar este análisis profundizando en la complejidad de la no linealidad, mejorando o ampliando el experimento MC y/o introduciendo un tercer parámetro que flexibilice el criterio RIC en la forma de t2>τ.log(q.r). De hecho, Foster y George (1994) también sugieren que el valor de τ óptimo puede diferir de 2 dependiendo justamente del tipo de colinealidad. Aplicando las fórmulas de inicialización de parámetros estimadas de acuerdo a lo previamente mencionado, se efectuó un extenso ejemplo, también mediante simulaciones MC, que permitió estudiar algunos aspectos sobre el desempeño del mismo. En general, puede afirmarse que el desempeño de algoritmo es bastante satisfactorio. Como es de esperar, mientras más baja es la colinealidad, más alto es el R2 ajustado

27

de la especificación verdadera, y/o mayor es el tamaño de muestra, mejores resultados se podrán esperar de este algoritmo. Las pruebas realizadas en este trabajo sugieren el tratar de asegurar un número suficientemente grande de iteraciones sin mejora (n). En este sentido, se propone que el investigador utilice como una convención informal a n ≥ 1500+10(q.r). APÉNDICE A Descripción de la subrutina GENDATA

I. Generación de los términos explicativos

(A.1) (A.2) (A.3) (A.4) (A.5) (A.6) Nota: para cada i =2,…,9 se generan diferentes conjuntos aleatorios según las expresiones desde (A.1) hasta la (A.3).

28

II. Generación de la variable dependiente

(A.7) (A.8) (A.9) (A.10) (A.11) (A.12)

29

APÉNDICE B Cuadro B.1

30

Cuadro B.2

31

Cuadro B.3

32

Cuadro B.4

33

REFERENCIAS BIBLIOGRÁFICAS

AKAIKE, H. (1973). “Information theory and an extension of the maximum likelihood principle”, in 2nd International Symposium on Information Theory, eds. B. N. Petrov and F. Csaki, Budapest: Akademia Kiado, pp. 267-281. BENJAMINI, Y., y HOCHBEERG, Y. (1995). “Controlling the false discovery rate: A practical and powerful approach to multiple testing”, Journal of the Royal Statistical Society, Ser. B, 57, 289-300. CLYDE, M., y GEORGE, E. I. (1999). “Empirical bayes estimation in wavelet nonparametric regression”, in Bayesian Inference in WaveletBased Models, eds. P. Muller and B. Vidakovic, New York: SpringerVerlag, pp. 309-322. CLYDE, M., y GEORGE, E. I. (2000). “Flexible empirical bayes estimation for wavelets”, Journal of the Royal Statistical Society, Ser. B, 681-698. DONOHO, D. L., y JOHNSTONE, I. M. (1994). “Ideal Spatial Adaptation by Wavelet Shrinkage,” Biometrika, 81, 425-256. FOSTER, D. P. y E. I. GEORGE (1994). “The Risk Inflation Criterion for Multiple Regression”, The Annals of Statistics, vol. 22, no 4, 19471975. GEORGE, E. I. y R. E. MCCULLOCH (1993). “Variable selection via Gibbs Sampling”, Journal of the American Statistical Association, vol. 88, no 423, 881-889. GEORGE, E. I. (2000). “The variable selection problem”, Journal of the American Statistical Association, vol. 95, no 452, 1304-1308.

34

GEORGE, E. I., y FOSTER, D. P. (2000). “Calibration and empirical bayes variable selection”, Biometrika, 87. HURVICH, C. M. y TSAI, C. L. (1989). “Regression and time series model selection in small samples”, Biometrika, 76, 297-307. MALLOWS, C. L. (1973). “Some comments on Cp”, Technometrics, 15, 661-676. RAO, C. R. y WU, Y. (1989). “A strongly consistent procedure for model selection in a regression problem”, Biometrika, 76, 369-374. SCHWARZ, G. (1978). “Estimating the dimension of a model”, The Annals of Statistics, 6, 461-464. SHAO, J. (1977). “An asymptotic theory for linear model selection”, Statistica Sinica, 7(2), 221-264. TIBSHIRANI, R., y KNIGHT, K. (1999). “The covariance inflation criterion for model selection”, Journal of the Royal Statistical Society, Ser. B, 61, 529-546. WEI, C. Z. (1992). “On predictive least squares principles,” The Annals of Statistics, 29, 1-42. YE, J. (1998). “On measuring and correcting the effects of data mining and model selection,” Journal of the American Statistical Association, 93, 120-131. ZHENG, X. y LOH, W. Y. (1997). “A consistent variable selection criterion for linear models with high-dimensional covariates”, Statistica Sinica, 7(2), 311-325.

35

Este N° 77 de la serie Documentos de Trabajo, en edición de 25 ejemplares, se terminó de imprimir en los Talleres de impresión del BCV, durante el mes de julio de dos mil ocho.