Modelo de regresión de Cox con métodos flexibles en pacientes con linfoma no hodgkin

Máster en Estadística e Investigación Operativa Título: Modelo de Regresión de Cox con Métodos Flexibles en Pacientes con Linfoma No Hodgkin Autor: Claudio Jaime Flores Flores Director: Guadalupe Gómez Melis Departamento: Departamento de Estadística e Investigación Operativa Convocatoria: 31 / enero / 2011 Universitat Polit`ecnica de Catalunya Facultat de Matem`atiques i Estad´ıstica Tesis de Master MODELO DE REGRESI ´ ON DE COX CON METODOS FLEXIBLE EN PACIENTES CON LINFOMA NO HODGKIN Claudio Jaime Flores Flores Director: Guadalupe G´omez Melis Departamento de Estad´ıstica e Investigaci´on Operativa A mis padres: Alejandro (en memor´ıa) y Br´ıgida Resumen En muchos estudios cl´ınicos es muy frecuente el uso de modelo de riesgo proporcional de Cox; el cual asume riesgos proporcionales y restringe a que el logaritmo de la raz´on de riesgo sea lineal en las covariables, lo cual en muchos casos no se verifica. En este sentido, una forma funcional no lineal del efecto de las covariables puede ser aproximada por una funci´on spline. En este trabajo, se presenta una aplicaci´on de los m´etodos flexibles como P-splines y polinomio fraccional para determinar y aproximar la forma funcional del efecto de las covariables (factores pron´osticos) en la supervivencia global de los pacientes con LNH. Los resultados muestran que el efecto de las covariables continuas como Hb, Leucocitos, linfocitos y la DHL presentan una forma funcional no lineal con el logaritmo de la raz´on de riesgo. Palabras clave: Modelo de Cox, P-spline, Polinomio Fraccional, LNH. MSC2000: Codis de la Mathematic Subject Classification Abstract In many clinical studies, Cox proportional hazard model is very common to use, it assumes proportional hazard and restricts the log hazard ratio to be linear in the covariates; these asumptions can not be verified. In this way, a nonlinear functional form of the covariates effect can be approximated by a spline function. In this paper, we present an application of flexible methods such as P-spline and fractional polynomial to identify, and align the functional form of the effect of covariates (prognostic factors) in the overall survival of patients with NHL. These results show that the effect of continuous covariates as: Hb, leukocytes, lymphocytes and LDH have a nonlinear functional form with the log hazard ratio. Keywords: Cox model, P-spline, polynomial fractional, NHL. MSC2000: 2000Mathematical Subject Classification ´ Indice general Cap´ıtulo 1. Introducci´on 1 1.1. Introducci´on 1 1.2. Objetivos del trabajo 2 1.3. Detalles del desarrollo 3 Cap´ıtulo 2. LNH: Aspectos cl´ınicos y metodol´ogicos 5 2.1. Linfoma No Hodgkin 5 2.2. Supervivencia y pron´ostico 6 2.3. Factores pron´osticos 7 2.4. Descripci´on de los datos 8 Cap´ıtulo 3. Conceptos b´asicos en an´alisis de supervivencia 11 3.1. Datos en an´alisis de supervivencia 11 3.2. Funciones del tiempo de supervivencia 14 3.3. Funci´on de m´axima verosimilitud 15 3.4. Procesos de conteo 16 3.5. Suavizamiento con splines 19 Cap´ıtulo 4. Modelo de regresi´on de Cox con m´etodos flexibles 27 4.1. Revisi´on de la literatura 27 4.2. Modelo de regresi´on de Cox 29 4.3. Modelo de regresi´on de Cox con P-splines 33 4.4. Modelo de regresi´on de Cox con polinomio fraccional 38 4.5. M´etodos de diagn´ostico en el modelo de Cox 39 Cap´ıtulo 5. Factores pron´osticos en LNH 47 5.1. Descripci´on de los datos 47 5.2. Aplicando el modelo de Cox cl´asico 48 5.3. Aplicando el modelo de Cox con P-splines 52 5.4. Aplicando el modelo de Cox con PF 56 5.5. Comparaci´on de los modelos 58 Cap´ıtulo 6. Discusi´on y conclusiones 63 Bibliograf´ıa 65 i 2.3. FACTORES PRON ´ OSTICOS 7 2.3. Factores pron´osticos Desde la d´ecada de los a˜nos 70 ha existido un enorme inter´es en la comunidad cient´ıfica por el estudio de los factores pron´osticos en los linfomas, como prototipo de enfermedad curable. Probablemente los linfomas sean las neoplasias mejor y m´as ampliamente estudiadas; sin embargo, se precisan de nuevos estudios para clarificar su utilidad real, debido a la aparici´on de nuevos FP (hemoglobina, leucocitos, linfocitos y marcadores tumorales), la existencia de un n´umero relativamente importante de pacientes que presentan recurrencia o que fallecen a consecuencia de la enfermedad y el desarrollo de nuevos m´etodos estad´ısticos para su an´alisis. Seg´un la literatura, los factores pron´osticos en los pacientes con LNH se agrupan en tres grandes grupos, aquellos que se derivan de las caracter´ısticas del paciente, del tumor y del tratamiento. Dentro de los FP dependientes del tumor, se tiene en cuenta las caracter´ısticas biol´ogicas del mismo y la carga tumoral (Mounier, et al. (1997), Costas, et al. (1998), Horsman y Hancock (2001), Rabasa 2002)). Dentro de los factores pron´osticos dependientes del paciente se considera la edad, estado funcional, enfermedades preexistentes y la competencia inmunol´ogica. La edad se considera como un factor pron´ostico, debido a que esta se asocia a una mayor morbimortalidad despu´es de los 60 ´o 70 a˜nos. La capacidad funcional del paciente, seg´un la escala ECOG (Zubrod), se considera como un factor pron´ostico al valorar la repercusi´on que la enfermedad produce en el estado general del paciente. Todas las situaciones cl´ınicas previas del paciente que puedan influir en la morbimortalidad y tolerancia al tratamiento se consideran como factores pron´osticos (ejemplo, las enfermedades cardiovasculares, diabetes, hepatitis, etc.). Los linfomas que aparecen en situaciones de inmunodeficiencia tienen un curso m´as agresivo y peor pron´ostico; prueba de ello son los LNH asociados al s´ındrome de inmunodeficiencia adquirida (SIDA). Dentro de los factores pron´osticos dependientes del tumor se considera el subtipo histol´ogico, el inmunofenotipo (c´elulas B o T), las alteraciones citogen´eticas, actividad proliferativa, extensi´on de la enfermedad y otras variables con significado pron´ostico. El subtipo histol´ogico, el patr´on de infiltraci´on (folicular o difuso) y aspecto citol´ogico de las c´elulas as´ı como su diferenciaci´on que dividen a los LNH en grupos diferenciados se consideran como factores pron´osticos; sin embargo, el pron´ostico de las diferentes entidades no son sustancialmente diferentes. El estudio del inmunofenotipo, ha permitido establecer el significado pron´ostico, al demostrarse un peor pron´ostico del linaje T frente al B. Las anomal´ıas citogen´eticas est´an presentes en la mayor´ıa de los LNH, la presencia de alteraciones cromos´omicas y el n´umero de ´estas reviste peor pron´ostico, mientras la ausencia se ha visto asociada a una mayor supervivencia. La extensi´on de la enfermedad, definida como la cantidad del tumor al momento del diagn´ostico, reviste una importancia pr´onostica en los LNH. Los siguientes par´ametros son considerados como factores pron´osticos: el estadio cl´ınico, n´umero y localizaci´on de ´areas ganglionares y extraganglionares afectas, tama˜no del tumor (”masa Bulky”, aquella masa cuyo tama˜no del di´ametro mayor es superior a 10cm), carga tumoral (n´umero de regiones ganglionares extensas (bulky) y el n´umero de localizaciones extraganglionares). 8 2. LNH: ASPECTOS CL´ INICOS Y METODOL ´ OGICOS Otras variables con significado pron´ostico son: presencia de s´ıntomas B (fiebre, sudoraci´on nocturna y p´erdida de peso), hemoglobina baja, leucocitos elevados, deshidrogenasa l´actica y β2-microglobulina elevada. 2.4. Descripci´on de los datos La base de datos objeto de an´alisis, corresponde a los datos de 2160 pacientes mayores o iguales a 14 a˜nos de edad con diagn´ostico de LNH, que fueron diagnosticados y tratados en el Instituto Nacional de Enfermedades Neopl´asicas (INEN), LimaPer´u, entre 1990 y 2002. As´ı mismo, cabe resaltar que los datos corresponden a una sub-base del estudio retrospectivo cl´ınico, patol´ogico y epidemiol´ogico de los LNH. El tratamiento que hab´ıan recibido los pacientes, seg´un la practica cl´ınica habitual, fue quimioterapia en la mayor´ıa de los casos (91.2 %) y en los restantes (8.8 %) radioterapia y/o cirug´ıa. El esquema de quimioterap´ıa que hab´ıan recibido fue generalmente (81.6 %) CHOP (ciclofosfamida, doxorubicina, vincristina y prednisona). 2.4.1. Descripci´on de los datos. Los datos recopilados al diagn´ostico de las variables relacionadas a las caracter´ısticas del paciente y del tumor (caracter´ısticas cl´ınicas) fueron: Edad: en a˜nos. G´enero: Femeno o masculino. Zubrod: Estado funcional del paciente, seg´un la escala ECOG. Primario: Localizaci´on ganglionar o extraganglionar. Tumor: Di´ametro mayor del tumor. Numero de ganglios afectados. Numero de sitios extraganglionares. Estadio cl´ınico: Extensi´on de la enfermedad, seg´un la clasificaci´on Ann Arbor. Sitios de met´astasis: Extensi´on del tumor a otras regiones o ´organos. S´ıntomas: Fiebre, sudoraci´on nocturna o baja de peso sin causa alguna. Tipo de LNH: Subtipo histol´ogico, seg´un la clasificaci´on disponible. VIH/SIDA: Infecci´on por VIH o SIDA. Hemoglobina: en g/dl. Leucocitos: N´umero de leucocitos /mm3. Linfocitos: Linfocitos en porcentaje. Deshidrogenasa l´actica: en UI/L. β2-microglobulina: en mg/L De los cuales, no se incluye en el an´alisis la informaci´on de las siguientes variables: tama˜no del tumor, n´umero de ganglios y el sitio de met´astasis, debido a que estas variables ya est´an reflejadas en el estadio cl´ınico. As´ı mismo, no se incluye el subtipo histol´ogico y el genotipo (c´elulas T o B) debido a que los pacientes fueron diagnosticados con tres criterios de clasificaci´on histopatol´ogica diferentes (Rappaport y Kiel, formulaci´on de trabajo (WF) y la clasificaci´on REAL). Y la β2M, debido a que este dato hab´ıa sido solicitado en muy pocos pacientes. 2.4. DESCRIPCI ´ ON DE LOS DATOS 9 En la Tabla 1 se muestra el n´umero de casos por cada una de las categor´ıas de las variables evaluables. El 36.9 % eran mayores de 60 a˜nos de edad, 50.9 % eran de sexo masculino, 27.0 % presentaban zubrod entre 2-4, 66.7 % ten´ıan enfermedad ganglionar, 49.2 % presentaban enfermedad en estadio cl´ınico avanzado (EC IIIIV), 38.1 % hab´ıan presentado s´ıntomas B, 48.1 % hab´ıan presentado hemoglobina bajo (Hb <12g/dl), 17.7 % leucocitos elevados (leucocitos >10mil), 7.7 % linfocitos elevados (linfocitos >40 %) y 60.0 % niveles de DHL elevados (>240 UI/L). Tabla 1. Variables y n´umero de casos por cada categor´ıa. Variables Mediana/rango Categor´ıa Casos Porcentaje ( %) Edad (a˜nos) 54.0 / 14 - 96 ≤60 1363 63.1 >60 797 36.9 G´enero - Femenino 1061 49.1 Masculino 1099 50.9 Zubrod - 0-1 1577 73.0 2-4 583 27.0 Primario - Extragang 719 33.3 Ganglionar 1441 66.7 Estadio cl´ınico - I-II 1097 50.8 III-IV 1063 49.2 S´ıntomas - A 1336 61.9 B 824 38.1 Hemoglobina (g/dl) 11.8 / 2.2 - 17.8 ≥12 1122 51.9 <12 1038 48.1 Leucocitos (103/mm3) 6.7 / 0.560 - 163 ≤1031777 82.3 >103383 17.7 Linfocitos ( %) 20 / 1 - 94 ≤40 1993 82.3 >40 167 7.7 Deshidrogenasa l´actica (UI/L) 298 / 24 - 8440 ≤240 863 40.0 >240 1297 60.0 Por otro lado, en la Gr´afica 1 se muestra la distribuci´on de las covariables continuas y se observa que todas ellas presentan asimetr´ıa. La asimetr´ıa es m´as pronunciada para las variables leucocitos y la deshidrogenasa l´actica (DHL); las cuales, son posibles debido a que los valores de los leucocitos pueden var´ıar entre <1000 y 100mil /mm3y de los DHL entre <240 y 10mil UI/L. 10 2. LNH: ASPECTOS CL´ INICOS Y METODOL ´ OGICOS edad Density 0 20 40 60 80 100 0.000 0.005 0.010 0.015 0.020 Hb Density 0 5 10 15 20 0.00 0.05 0.10 0.15 0.20 leucocitos Density 0 20000 40000 60000 80000 0.00000 0.00005 0.00010 0.00015 linfócitos Density 0 20 40 60 80 100 0.000 0.010 0.020 0.030 DHL Density 0 2000 6000 10000 0.0000 0.0005 0.0010 0.0015 0.0020 Gr´ afica 1. Distribuci´on de las variables continuas de los datos de pacientes con LNH. Cap´ıtulo 3 Conceptos b´asicos en an´alisis de supervivencia En este cap´ıtulo se presentan algunos conceptos b´asicos en an´alisis de supervivencia, que ser´an utilizadas en los siguientes cap´ıtulos como son: datos de supervivencia, funciones de distribuci´on, procesos de conteo y el m´etodo splines. 3.1. Datos en an´alisis de supervivencia En muchos estudios cl´ınicos los investigadores pueden estar interesados en analizar el tiempo de supervivencia y las caracter´ısticas cl´ınicas de los pacientes que pueden estar relacionados con el tiempo de supervivencia. En este contexto, hay tres aspectos caracter´ısticos en el an´alisis de los datos de supervivencia, que son: El tiempo de seguimiento hasta la ocurrencia del evento, llamado tiempo de supervivencia. La censura o censuramiento, que se origina debido a estudios que son terminados antes de los resultados de todas las unidades conocidas, y La presencia de las variables explicativas, llamadas covariables. 3.1.1. Tiempo de supervivencia. Se denomina tiempo de supervivencia al tiempo transcurrido entre la fecha de inicio (ingreso al estudio, inicio de tratamiento, etc.) y la fecha de ocurrencia de un evento de inter´es (reca´ıda, progresi´on, muerte, etc.) o la fecha o tiempo en que finaliza el estudio. En general, el tiempo de supervivencia es un proceso continuo, la longitud de la duraci´on puede medirse utilizando un n´umero real no negativo. Como un t´ermino gen´erico, el tiempo desde la iniciaci´on de un evento (nacimiento, diagn´ostico, inicio de tratamiento, etc) hasta la ocurrencia de un evento de inter´es (reca´ıda, progesi´on, muerte, etc.) se denomina como tiempo de supervivencia, a´un cuando el evento final es algo diferente de la muerte. Algunos ejemplos: Tiempo desde el nacimiento hasta la muerte. Tiempo desde el nacimiento hasta el diagn´ostico de c´ancer. 11 12 3. CONCEPTOS B ´ ASICOS EN AN ´ ALISIS DE SUPERVIVENCIA Tiempo transcurrido desde la aparici´on de la enfermedad hasta la muerte. Tiempo transcurrido desde la respuesta cl´ınica a la reca´ıda. En general, para fines de nuestra aplicaci´on, se considera el an´alisis de supervivencia cl´asico que se centra en el tiempo hasta la ocurrencia de un simple evento (muerte del paciente) para cada individuo, o m´as exactamente el tiempo transcurrido desde inicio hasta la ocurrencia del evento muerte. 3.1.2. Censura o Censuramiento. Normalmente, los estudios de supervivencia tienen una duraci´on predeterminada, por lo que no todos los sujetos en seguimiento habr´an fallado a su finalizaci´on. Por lo tanto, el investigador sabr´a que un cierto n´umero de individuos han ”sobrevivido”durante el periodo de tiempo, pero desconocer´a el momento exacto en que hubiera fallado si el estudio si hubiera prolongado de forma indefinida. A este tipo de datos se llaman observaciones censuradas. Se dice que las observaciones est´an censuradas cuando contienen informaci´on parcial sobre los tiempos de supervivenvia durante un periodo de seguimiento. La informaci´on parcial, ocurre debido a causas como: Retiro del estudio por causas ajenos al evento de inter´es P´erdida de acompa˜namiento, o Finalizaci´on del estudio. En general, el t´ermino de censura hace referencia a un tipo de p´erdida de informaci´on en situaciones en las que la variable de inter´es es el tiempo de supervivencia. La censura surge en las ocasiones en las que hay individuos de la muestra para los que no se conoce exactamente su tiempo de supervivencia, sino que ´unicamente se sabe que ´este ha ocurrido dentro de un cierto intervalo de tiempo. De esta forma se puede considerar tres tipos de censura: censura por la derecha, por la izquierda y censura en un intervalo. Se dice censura por la derecha, cuando en el momento en que finaliza el estudio hay sujetos para los que no se conoce el instante exacto de falla, sino que solamente se sabe que es posterior a un momento dado. Censura por la izquierda, cuando el momento exacto en que ocurri´o la falla es desconocido, tan s´olo se sabe que ha ocurrido antes de que el sujeto se incluya en el estudio. Y censura en un intervalo cuando el evento de inter´es no se puede observar exactamente y s´olo se sabe que ha ocurrido en un cierto intervalo de tiempo. En los casos de censura por la derecha se tiene tres tipos de censura: censura de tipo I (censuramiento por tiempo), tipo II (censuramiento por fallas) y tipo III (censuramiento aleatorio). Si un estudio termina en un tiempo pre-establecido y algunos de los tiempos de supervivenvia son no observados, tenemos censura tipo I. En el caso que el estudio termina despu´es de la ocurrencia de una determinada cantidad pre-establacida de eventos, tenemos censura de tipo II. Para prop´ositos del presente trabajo nos enfocaremos en datos con censura o censurados por la derecha y de tipo aleatorio. La censura aleatoria surge de manera 3.1. DATOS EN AN ´ ALISIS DE SUPERVIVENCIA 13 natural en las investigaciones biom´edicas debido a que los pacientes entran al estudio en tiempos diferentes y de manera aleatoria, y que cada paciente tiene un modo propio de censura, debido a cualquiera de las causas descritas (retiro, p´erdida de siguimiento, finalizaci´on del estudio), de modo que los tiempos de censuramiento son tambi´en aleatorias. 3.1.3. Covariables. Adem´as de los datos de supervivencia (tiempo de supervivencia y la variable indicadora de censura), tambi´en se pueden observar otros datos, variables que representan la heterogeneidad existente en la poblaci´on, tales como, la edad, el g´enero, la hemoglobina, estadio cl´ınico, etc. Estas variables son conocidas como variables explicativas o covariables y son muy frecuente en muchos estudios cl´ınicos. En general, las covariables son variables independientes y observables. Seg´un la escala de medici´on, se pueden clasificar en variables cuantivativas y cualitativas, y segun la evoluci´on en el tiempo en variables fijas o tiempo-dependientes. 3.1.3.1. Clasificac´ıon seg´un la escala de medici´on. Las variables cuantitativas son variables que se pueden medir expres´andose num´ericamente. Estas variables pueden ser de dos tipos: Continuas, cuando admiten tomar cualquier valor dentro de un rango num´erico (ejemplo: edad, peso, talla, tama˜no del tumor, hemoglobina, etc). Discretas, si solamente toman valores enteros, por lo que no admiten los valores intermedios en un rango dado (ejemplo: n´umero de hijos, n´umero de ganglios, etc). Las variables cualitativas son variables que representan distintas cualidades de un individuo; estas cualidades se denominan atributos o categor´ıas, y la medici´on de ´estas variables consiste en la clasificaci´on de dichos atributos. En el proceso de medici´on de las variables cualitativas, se pueden utilizar dos escalas, la ordinal y la nominal. En la escala ordinal, la clasificaci´on de las categor´ıas presentan un orden natural (ejemplo: grupos de edad, estado de performance, estadio cl´ınico, etc). En la escala nominal, las categor´ıas de la variables no se clasifican de acuerdo a un criterio de orden tanto inherente como jer´arquico (ejemplo: g´enero, estado civil, etc). Dependiendo de los valores que tome una variable cualitativa, ´estas pueden ser dicot´omicas o bien polit´omicas. 3.1.3.2. Clasificac´ıon seg´un la evoluci´on en el tiempo. Las variables fijas, son variables cuyos valores no var´ıan durante la evolucion del estudio; es decir el valor de estas variables no cambian durante el periodo de seguimiento. El valor o el atributo de estas variables al inicio del estudio es la misma en cualquier momento del tiempo (ejemplo: g´enero, raza, tipo de Rh, etc). Una complicaci´on que puede ocurrir en el an´alisis de supervivencia es observar variables que pueden variar en el tiempo, denominadas variables tiempo-dependientes. Los valores de estas variables no es lo mismo al final que al inicio del estudio. (ejemplo: edad, estado de la enfermedad, respuesta al tratamiento, modificaci´on de la dosis de un determinado medicamento a lo largo del tratamiento, etc.) 14 3. CONCEPTOS B ´ ASICOS EN AN ´ ALISIS DE SUPERVIVENCIA 3.2. Funciones del tiempo de supervivencia En an´alisis de supervivencia existen dos funciones de gran inter´es: la funci´on de supervivencia y la funci´on de riesgo; las cuales, son descritas brevemente en esta secci´on. Para comenzar con el modelado de datos de supervivencia se partir´a de la suposici´on de que la poblaci´on es homog´enea. Sea Tel tiempo de supervivencia, una variable aleatoria positiva (T > 0) con funci´on de densidad f(t) y funci´on de distribuci´on F(t). La funci´on de densidad es definida como el l´ımite de la probabilidad de un individuo de fallecer en un intervalo de tiempo [t, t +4t) por unidad de tiempo, y es expresada por f(t) = l´ım 4t→0 P(t≤T < 4t) 4t(1) y la funci´on de distribuci´on es definida como la probabilidad de fallecer de un invididuo antes de un tiempo t, y es expresada por F(t) = P(T≤t) Adem´as de fyF, en el an´alisis de supervivencia, la distribuci´on del tiempo de supervivencia puede ser caracterizada por otras funciones equivalentes: la funci´on de supervivencia y la funci´on de riesgo. La funci´on de supervivencia es definida como la probabilidad de que un individuo sobreviva por lo menos un determinado tiempo t. Esta funci´on es decreciente con un valor 1 para T= 0 y cero para T=∞, y es expresada como S(t) = P(T > t) = 1 −F(t) = 1 −Z∞ t f(s)ds. (2) La funci´on de riesgo se define como la probabilidad condicional de que un sujeto muera en un intervalo de tiempo (t, t + ∆t) dado que ya sobrevivi´o por lo menos un tiempo t, interpretado como la tasa instant´anea de falla, y es expresada como λ(t) = l´ım ∆t→0 P(t≤T < t + ∆t/T > t) ∆t De la expresi´on (1) y (2) la funci´on de riesgo es, λ(t) = f(t)/S(t). Por otro lado, la funci´on de riesgo acumulado se denota como H(t) = R∞ tλ(s)ds, por lo tanto, la funci´on de supervivencia puede ser calculada a partir de la funci´on de riesgo por S(t) = exp(−H(t)) = exp −Z∞ t λ(s)ds En consecuencia, en el an´alisis de los datos de supervivencia, la descripci´on y modelado de los tiempos de supervivencia se puede realizar con cualquiera de estas funciones. 3.3. FUNCI ´ ON DE M ´ AXIMA VEROSIMILITUD 15 3.3. Funci´on de m´axima verosimilitud Una de las partes fundamentales de todo procedimiento estad´ıstico, es la estimaci´on de los par´ametros del modelo estad´ıstico, basado en una muestra. El procedimiento de estimaci´on en una muestra no censurada no es complicado; sin embargo, en una muestra censurada el procedimiento de estimaci´on de los par´ametros esta sujeto a un factor o indicador de censura o censuramiento. En las investigaciones biom´edicas, el censuramiento surge de manera natural, debido a que los pacientes entran al estudio en tiempos diferentes, de manera que cada paciente tiene un modo propio de censuramiento, debido a cualquiera de las tres causas descritas anteriormente (perdida de acompa˜namiento, retiro del estudio por eventos ajenos al estudio y t´ermino del estudio), de modo que los tiempos de censuramiento son tambi´en aleatorios. En esta secci´on se describe brevemente la funci´on de verosimilitud para el procedimiento de estimaci´on, mediante el m´etodo de m´axima verosimilitud para los par´ametros del modelo o la funci´on de supervivencia bajo censuramiento por la derecha y de manera aleatoria. Sea Yel tiempo de supervivencia y Cel tiempo de censuramiento asociado, con funci´on de densidad y funci´on de supervivencia, (fT(t), ST(t)) y (gC(t) y 1−GC(t)), respectivamente. Bajo un mecanismo de censura no-informativo, se supone que Y yCson independientes. Tambi´en se asume que G(t) no depende de ninguno de los par´ametros de S(t), por lo que no aporta informaci´on alguna para la distribuci´on del tiempo de supervivencia. En este modelo de censura, lo que se observa por unidad muestral es el par aleatorio (T, δ) definido como T= min(Y, C), y δ=I[Y≤C]=(1,si Tes no censurado 0,si Tes censurado donde δes la variable indicadora de censura. Sean los tiempos de supervivencia observados para nindividuos que consiste de los pares (t1, δ1),...,(tn, δn). La funci´on de verosimilitud es dada por L= n Y i=1 [f(ti)(1 −Gi)]δi[giS(ti)]1−δi.(3) Debido a que el tiempo de censuramiento es no informativo, la funci´on de verosimilitud se reduce en t´erminos de f(t) y S(t) a: L= n Y i=1 [f(ti)]δi[S(ti)]1−δi.(4) 16 3. CONCEPTOS B ´ ASICOS EN AN ´ ALISIS DE SUPERVIVENCIA En consecuencia, reemplazando las funciones respectivas en las expresiones (3) y (4), se pueden obtener los estimadores de los par´ametros del modelo de distribuci´on supuesto para variable aleatoria Ty las funciones equivalentes. 3.4. Procesos de conteo En el an´alisis de supervivencia, el enfoque del an´alisis est´a en la observaci´on de la ocurrencia de eventos sobre el tiempo. Dichas ocurrencias constituyen procesos puntuales. Estos procesos pueden ser descritos como el conteo del n´umero de eventos que se van presentando durante el tiempo, lo que lleva al t´ermino de ” procesos de conteo”. Algunos ejemplo pueden ser: Contar el n´umero de veces que una persona se despierta durante la noche. Contar las muertes en un grupo de pacientes con tratamiento en un ensayo cl´ınico. En general, hay una teor´ıa matem´atica rigurosa para los procesos de conteo, los cuales son extremadamente ´utiles para el an´alisis estad´ıstico de los datos de supervivencia y datos de eventos hist´oricos. La raz´on de usar procesos de conteo y martingalas es porque ´estos m´etodos proporcionan formas directas de estudiar las propiedades de muestras grandes de los estimadores. En esta secci´on se describe brevemente los t´erminos b´asicos del proceso de conteo. 3.4.1. Procesos de conteo. Los tiempos de supervivencia pueden ser representados a trav´es de ciertos procesos estoc´asticos. Los datos en s´ı pueden ser descritos como un proceso de conteo, el cual es simplemente una funci´on aleatoria del tiempo t,N(t). Esta funci´on es cero en el tiempo inicial y constante en el tiempo, excepto en el tiempo donde ocurre el evento, donde hace un salto de tama˜no 1. Considerando para un solo tipo de evento. Para un tiempo tdado, sea N(t) el n´umero de eventos que han ocurrido hasta el tiempo. Entonces N(t) es un proceso de conteo. Sean los datos de 10 pacientes de un estudio cl´ınico hipot´etico, escrito en orden creciente: 2.70, 3.50+, 3.80, 4.19, 4.42, 5.43, 6.32+, 6.46+, 7.32, 8.11+. El proceso de conteo correspondiente para estos datos se ilustra en la Gr´afica 2. En la Gr´afica 2 se observa que el proceso da saltos de una unidad en cada evento observado en el tiempo, es constante entre los eventos y es continua por la derecha. As´ı, los tiempos de supervivencia pueden ser representados a trav´es de ciertos procesos estoc´asticos. En los p´arrafos siguientes se describen los t´erminos b´asicos de uso com´un en la metodolog´ıa de los procesos de conteo. 3.5. SUAVIZAMIENTO CON SPLINES 23 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 x base[, i] Bases B−spline de grado 1 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 x base[, i] Bases B−spline de grado 3 Gr´ afica 3. Bases B-splines de grado 1 y 3 definidos en [0,1] con 9 nodos equi-espaciados, respectivamente. 3.5.4. P-splines. Eilers y Marx (1996) introducen el t´ermino P-splines, llamado splines penalizado (Ruppert y Carroll, 2000) o pseudo-splines (Hastie, 1996), son una extensi´on de B-splines y comparten muchas de las propiedades. Los P-splines es un t´ermino intermedio entre el suavizamiento splines y regresi´on splines, de hecho combinan lo mejor de ambos enfoques. Los P-splines utilizan menos par´ametros que los splines de suavizamiento, pero la selecci´on de los nodos no es tan determinante como en los splines de regresi´on. Son splines de rango bajo, el n´umero de nodos es mucho menor que la dimensi´on de los datos, al contrario de lo que ocurre en el caso de los splines de suavizamiento. El n´umero de nodos, en el caso de los P-splines, no supera los 40, lo que hace que sean computacionalmente m´as eficientes, sobre todo cuando se trabaja con gran cantidad de datos. Adem´as, la introducci´on de penalizaciones relaja la importancia de la elecci´on del n´umero y la localizaci´on de los nodos, cuesti´on que es de gran importancia en los splines de rango bajo sin penalizaciones (Rice y Wu, 2001). La novedad es que los autores proponen el uso de B-splines sim´etrico y penalizan estos, no en la segunda derivada, sino en las diferencias entre los coeficientes splines adyacentes. Este tipo de penalizaci´on es m´as flexible ya que es independiente del grado del polinomio utilizado para construir los B-splines, un criterio que es f´acil de implementar, que resulta estrechamente relacionado con la usual penalizaci´on. Sean los t´erminos flexibles expresados como: 24 3. CONCEPTOS B ´ ASICOS EN AN ´ ALISIS DE SUPERVIVENCIA f(x) = M X m=1 αmBm(x),(9) donde Bm(x) son las funciones bases B-splines. El tipo de penalizaci´on est´a estrechamente ligado al tipo de base que se utilice. Si se utilizan polinomios truncados, la penalizaci´on es ”ridge” independientemente del grado de los polinomio truncados. Para las bases B-splines el t´ermino de penalizaci´on frecuentemente utilizada es la integral de la segunda derivada al cuadrado de la curva (funci´on B-spline), tal como fue sugerida por O’Sullivan (1986), (1/2)λZ[f00(x)]2dx. (10) Supongamos que tenemos una base B-splines construida con knodos y utilizamos m´ınimos cuadrados penalizados para ajustar el modelo. La Gr´afica 4 muestra ajuste de una curva mediante B-splines sin y con penalizaci´on, las funciones que forman la base multiplicadas por los coeficientes, as´ı como los coeficientes (en c´ırculo). El efecto de la penalizaci´on que fuerza a los coeficientes a seguir un patr´on m´as suave. En la parte izquierda se aprecia como el patr´on err´atico de los coeficientes da lugar a una curva poco suave, en cambio en la parte derecha, cuando se impone que se pase de un coeficiente a otro de forma suave, la curva tambi´en lo es. 3.5. SUAVIZAMIENTO CON SPLINES 25 0.0 0.2 0.4 0.6 0.8 1.0 −1 0 1 2 x y 0.0 0.2 0.4 0.6 0.8 1.0 −1 0 1 2 x y Gr´ afica 4. Curva estimada con 20 nodos mediante las bases B-splines sin y con penalizaci´on. Cap´ıtulo 4 Modelo de regresi´on de Cox con m´etodos flexibles En este cap´ıtulo se presenta una breve revisi´on de literatura y se describe la base te´orica del modelo de regresi´on de Cox utilizando los m´etodos flexibles mediante los splines penalizados y el polinomio fraccional. As´ı mismo, se describe los procedimientos de diagn´ostico del modelo de Cox. 4.1. Revisi´on de la literatura En muchos trabajos de investigaci´on m´edica es muy com´un que en cada paciente adem´as de ser observado el tiempo de supervivencia sean observadas las caracter´ısticas cl´ınicas de los pacientes, llamadas covariables. Si el inter´es es determinar el efecto de las covariables en el tiempo de supervivencia, el prop´osito del estudio se centra en el an´alisis de las relaciones entre el tiempo de supervivencia y las covariables mediante un modelo de regresi´on. Los modelos parametricos son eficientes cuando se tiene informaci´on del modelo de distribuci´on subyacente a las variables y s´olo resta por determinar un n´umero finito de par´ametros; sin embargo, una fuente de error puede ser elegir una familia param´etrica no adecuada. En estos casos podemos utilizar los modelos semiparam´etricos o no-param´etricos que adem´as de permitir graduar las probabilidades brutas que no siguen un modelo param´etrico establecido, pueden utilizarse para proporcionar una prueba de diagn´ostico de los modelos param´etricos o simplemente para explorar los datos. En an´alisis de datos de supervivencia, el modelo de regresi´on semiparam´etrica muy frecuentemente utilizado es el modelo de riesgos proporcionales de Cox (Cox, 1972), llamado generalmente como modelo de Cox. Este modelo asume que la funci´on de riesgo es constante sobre un periodo de tiempo y el efecto de las covariables se relaciona linealmente con el logaritmo de la raz´on de riesgos. Si los supuestos del modelo no se cumplen, el modelo de Cox no es el m´as adecuado, entonces el modelo de Cox estratificado o el modelo de Cox con variables tiempo-dependiente podr´ıan ser una alternativa. Otras alternativas pueden ser el modelo de odds proporcional y el modelo log-log´ıstico. 27 28 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES Sin embargo, el modelo de Cox estratificado no es adecuado cuando si se tiene m´as de una variable que no verifica el supuesto de riesgos proporcionales, ya que se pierde la informaci´on para estas variables que pueden ser de importancia para qui´en investiga. En el modelo de Cox con variable tiempo-dependiente, as´ı como en los dem´as modelos alternativos si el efecto de las covariables presentar´a una forma funcional no-lineal, el problema de modelado a´un estar´ıa pendiente. Durante las dos ´ultimas d´ecadas numerosas t´ecnicas se han desarrollado para explorar la forma funcional del efecto de las covariables de una manera m´as flexible, utilizando para ello m´etodos de suavizaci´on (smoothing spline) y polinomios fraccionales (fractional polynomials o FP). En este trabajo se utilizan estos dos m´etodos de modelamiento como una alternativa en situaciones en que la forma funcional del efecto de las covariables en la funci´on de riesgo es no lineal. Los m´etodos splines son una herramienta usual en muchos contextos estad´ısticos, permiten el manejo de relaciones no lineales complejas, dif´ıciles de calcular con modelos param´etricos convencionales. Los splines m´as utilizados son los splines penalizados (penalized splines o P-splines). Los P-splines (Eilers y Marx, 1996) aproximan una funci´on desconocida por un spline polinomial que puede ser escrito como una combinaci´on lineal de funciones bases splines. En cambio, los FP (Royston y Sauerbrei, 2008) aproximan la funci´on desconocida por una suma de trasformaciones potencia de las covariables, que son m´as flexibles que los polinomios ordinarios, como tal admiten potencias negativas y no enteras. Las diferentes aproximaciones mediante splines son bastante amplias, abarcando desde t´ecnicas de suavizaci´on splines (Hastie y Tibshirani (1990), Wahba (1990), Green y Silverman (1994)) con nodos fijos, hasta el uso de regresi´on splines con selecci´on adaptativa de los nodos (Friedman, 1999). Eilers y Marx (1996) propusieron el uso de splines penalizado (penalized splines, P-splines), un enfoque diferente que puede ser visto como un compromiso entre suavizamiento spline y regresi´on spline. En el modelo de Cox, O’Sullivan (1988) utiliza splines de suavizado para estimar el efecto no-lineal de la covariable. Sleeper y Harrington (1990) utilizan regresi´on splines con un n´umero reducido de nodos, donde los coeficientes son estimados utilizando el m´etodo est´andar sin funciones de penalizaci´on; sin embargo, son m´as sensibles al n´umero y localizaci´on de los nodos y por tanto, son m´as inestables. Gray (1992) utiliza splines penalizados (penalizad splines) en el modelo de Cox. Para los prop´ositos de este trabajo, aqu´ı se presentan las bases metodol´ogicas de los modelos de regresi´on de Cox con m´etodos flexibles para determinar los factores pron´osticos en la supervivencia global de los pacientes con LNH y se estima la forma funcional del efecto de las covariables en la raz´on de riesgo (hazard ratio). Previamente se hace una breve descripci´on del modelo de Cox (modelo de Cox cl´asico), luego se describe el modelo de Cox con P-spline y modelo de Cox con polinomio fraccional, y finalmente se realizan las comparaciones entre los dos m´etodos (P-splines y FPs), en describir la forma funcional no lineal de los efectos de las covariables en la funci´on de riesgo en lugar de simplemente determinar el efecto significativo de las covariables. 4.2. MODELO DE REGRESI ´ ON DE COX 29 4.2. Modelo de regresi´on de Cox Sean los datos observados en la muestra de la forma (T,δ,X), donde, Tes el tiempo de supervivencia observada, δel indicador de censura y Xlas covariables. Sea el objetivo del an´alisis evaluar el efecto de las covariables en el tiempo de supervivencia hasta la ocurrencia de un evento de inter´es (ejemplo: falla, muerte, recurrencia u otros). En esta situaci´on el modelo de regresi´on de Cox es una alternativa para analizar el efecto de las covariables en el tiempo de supervivencia. 4.2.1. Modelo de Cox cl´asico. El modelo de riesgos proporcionales introducido por Cox (1972), llamado modelo de regresi´on de Cox o modelo de Cox, es de la forma λ(t/X) = λ0(t) exp{X0β},(11) donde λ0(t) es la funci´on de riesgo basal cuando X= 0 (cuya distribuci´on es no especificada), β= (β1,...,βp) es un vector de par´ametros del modelo y X= (X1,...,Xp) son los vectores de covariables. En el modelo de Cox, la funci´on de riesgo es el producto de la funci´on de riesgo basal λ0(t) y un escalar exp{X0β}que s´olo depende de los par´ametros y las covariables. Las cuales, imponen las siguientes restricciones al modelo (11): la raz´on entre las funciones de riesgo para dos individuos es proporcional (riesgo proporcional) el logaritmo de la raz´on de riesgo es independiente del tiempo (riesgo constante) el logaritmo de la raz´on de riesgo se relaciona linealmente con las covariables (forma funcional lineal). En el modelo de la forma (11), el problema de modelamiento de los datos de supervivencia se reduce a la estimaci´on de los par´ametros βyλ0(t) condicionada a ˆ β, y en realizar la prueba de hip´otesis sobre los par´ametros, H0:β=β0, para evaluar el efecto de las covariables en la funci´on de riesgo. 4.2.2. Estimaci´on de los par´ametros. Sea una muestra de tama˜no n, donde los datos observados en la muestra son las ternas (Ti, δi, Xi). Asumimos que los tiempos de supervivencia son censurados por la derecha bajo un mecanismo de censura no-informativa. En el modelo de la forma (11) la estimaci´on de los par´ametros βdeben ser estimados a partir de las observaciones muestrales. La presencia de la componente no-param´etrica invalida el uso del m´etodo de m´axima verosimilitud. Para estimar el vector de par´ametros β, Cox (1972) propone el m´etodo de verosimilitud parcial. Supongamos que los datos est´an compuestos por nindividuos, de los cuales existen r tiempos de muertes diferentes, los restantes n−rtiempos se consideran censurados. 30 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES As´ı mismo, se supondr´a que solo un individuo muere en cada tiempo, es decir no hay empates. Sea t(1) < . . . < t(r)los distintos tiempos de muerte ordenados (r≤n), R(t(j)) = Rj el conjunto de individuos que se encuentran a riesgo en el tiempo t(j), es decir, conjunto de individuos que se encuentran con vida y sin censura antes de t(j), y sea X(j)el vector de covariables asociadas. Cox (1972) define que la funci´on de verosimilitud parcial para el modelo de riesgos proporcionales como L(β) = r Y j=1 exp(X0 (j)β) P l∈Rj exp(X0 lβ)(12) Equivalentemente incluy´endose toda la muestra, sea t(1) < . . . < t(n)las distintas observaciones de tiempos ordenados, δ(i)indicador de censura respectiva y X(i)el vector de covariables asociadas. Sea m(i)el evento (muerte) de un individuo ien el instante t(i)y sea R(t(i)) = R(i) el conjunto de individuos en riesgo en el tiempo t(i). Dada la funci´on de riesgo de la forma (11) para cada sujeto iy definida la probabilidad de observar una falla en un individuo ide la forma P(m(i)/t(i)R(t(i))) =   exp(X0 iβ).X l∈R(t(i)) exp(X0 lβ)   δ(i) , el logaritmo de la funci´on de verosimilitud parcial es dado por `(β) = ln[L(β)] = n X l=1 δ(i)  X0 (l)β−ln X l∈R(ti) exp(X0 lβ)  (13) Cuando los datos contienen tiempos observados empatados, la verosimilitud parcial (12) tiene que ser modificada de alguna forma. Se han propuesto varias aproximaciones para la funci´on de verosimilitud parcial en esta situaci´on, por ejemplo, Breslow (1974), Efron (1977) y Cox (1972). Sea t(1) < . . . < t(r)los rdistintos tiempos observados y ordenados. Sea djel n´umero de fallos observadas en t(j)yDj≡D(t(j)) = j1,...,jdjel conjunto de etiquetas de los individuos que fallan en tj. Sea sj=Pl∈DjXlyRjel conjunto de sub´ındices de individuos que se encuentran en riesgo antes de tj. La aproximaci´on de la verosimilitud parcial sugerida por Breslow (1974) considera que las djfallas al tiempo t(j)son distintos y ocurren secuencialmente. Cuando se tienen pocos empates, esta proporciona una muy buena aproximaci´on de la funci´on de verosimilitud parcial. La verosimilitud debida a Breslow (1974), en el caso de empates, es 4.2. MODELO DE REGRESI ´ ON DE COX 31 r Y j=1 exp(s0 jβ) "P l∈Rj exp(s0 lβ)#dj(14) Una aproximaci´on sugerida por Efron (1977) es r Y j=1 exp(S0 jβ) Qdj k=1 hPl∈Rjexp(S0 lβ)−(k−1)d−1 kPl∈Djexp(X0 lβ)i(15) Otra aproximaci´on sugerida por Cox (1972) es r Y j=1 exp(S0 (j)β) P l∈R(t(j);dj) exp(S0 lβ)(16) donde R(t(j);dj) denota el conjunto de todos los subconjuntos de djindividuos seleccionados del conjunto de riesgo R(tj) sin reemplazo. De este modo, si `∈ R(t(j);dj), ´este es de la forma `1,...,`dj. La verosimilitud parcial anterior es computacionalmente dif´ıcil si el n´umero de empates es grande. A partir de cualquiera de las expresi´on (12), (13), (14), (15) y (16), se pueden obtener los estimadores de m´axima verosimilitud parcial de los par´ametros β= (β1,...,βp) y se pueden realizar las pruebas de hip´otesis basadas en la distribuci´on asint´otica de los estimadores. El estimador de m´axima verosimilitud del vector de par´ametros se obtiene como una soluci´on del sistema de pecuaciones no lineales generadas por las derivadas parciales de `(β) respecto a los par´ametros, βj, d`(β)/dβj= 0 (j= 1,...,p) y utilizando procedimientos iterativos como el m´etodo de Newton-Raphson. 4.2.3. Prueba de hip´otesis e inferencia. Bajo las condiciones de regularidad (consistencia, distribuci´on asint´otica normal, eficiencia asint´otica), se dice que los estimadores de m´axima verosimilitud tienen distribuci´on asint´oticamente normal con media β= (β1,...,βp) y matriz de varianza y covarianza I∗(β). Dada las complicaciones en el c´alculo de I∗(β) es com´un utilizar en su reemplazo la matriz de informaci´on observada I(β) . Sea la funci´on de puntaje definida como la derivada parcial del logaritmo de la funci´on de verosimilitud con respecto a los par´ametros, Uj(β) = d`(β)/dβj(j= 1,...,p) y la matriz de informaci´on como el negativo de la derivada del puntaje de eficiencia con respecto a los par´ametros, I(β) = [−dU(β)/dβk] (k= 1,...,p). En el modelo de la forma (11) la prueba de hip´otesis global, H0:β=β0vs. H1: β6=β0, para una muestra de tama˜no nsuficientemente grande, bajo la distribuci´on asint´otica de los estimadores de m´axima verosimilitud parcial, se pueden realizar mediante tres estad´ısticos de prueba diferentes: la prueba de Wald, la raz´on de 32 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES verosimilitud y la prueba del score, que bajo la hip´otesis nula todas ellas tienen distribuci´on χ2con pgrados de libertad. Sean b= (b1,...,bp) los estimadores de m´axima verosimilitud del vector de par´ametros desconocidos, β= (β1,...,βp) y I(β) la matriz de informaci´on. Los estad´ısticos de prueba para contrastar la hip´otesis, H0:β=β0, basados en la distribuci´on asint´otica de los estimadores son: La prueba basada en la normalidad asint´otica de los estimadores, referida como la prueba Wald: χ2 W= (b−β0)0(I(b))(b−β0)∼χ2 p La prueba de raz´on de verosimilitud: χ2 LR = 2[`(b)−`(β0)] ∼χ2 p La prueba de score: esta prueba es basada en el puntaje de eficiencia U(β) = (U1(β),...,U1(β)). Para muestra grande, U(β) es asint´oticamente normal pvariada con media 0 y matriz de varianza y covarianza I(β), χ2 Sc =U(β0)0 I−1(β0)U(β0)∼χ2 p Por otro lado, debido a que los m´etodos de diagn´ostico del modelo est´an desarrollados en el contexto de proceso de conteo y teor´ıa de martingalas, en los p´arrafos siguientes se describe el modelo de Cox en el contexto de proceso de conteo. 4.2.4. Modelo de Cox en el contexto de procesos de conteo. El tratamiento de datos de supervivencia mediante procesos de conteo tiene sus or´ıgenes en el trabajo de Aalen (1978). Posteriormente Andersen y Gill (1982) integraron el modelo de Cox en el marco de procesos de conteo, generalizando de esta forma el tratamiento habitual de los modelos de supervivencia. Andersen y Gill (1982) extienden el modelo Cox en el contexto de procesos de conteo y obtienen pruebas martingala para las propiedades asint´oticas de los estimadores asociados (Martinussen and Scheike, 2006). Supongamos que observamos nobservaciones de la forma (Ti, δi, Xi), donde Ties el tiempo de supervivencia censurado por la derecha, δien indicador de censura, Xiel vector de covariables. El modelo de Cox, asume que la funci´on de intensidad es de la forma λ(t) = Y(t)λ0(t) exp(X0β),(17) donde Y(t) es una indicador de riesgo que es uno si el evento no ha ocurrido, λ0(t) es la funci´on de riesgo base no-param´etrico localmente integrable, X= (X1,...,Xp) es un vector de pcovariables y βes el vector de par´ametros del modelo. Sea una muestra de tama˜no nconteniendo datos de la forma (Ni(t), Yi(t), Xi)i= 1,...,n que son observados en un intervalo de tiempo [0, t], t < ∞, y que cada Ni(t) tiene intensidad de la forma (17). Los par´ametros βdel modelo son estimados maximizando como en la funci´on de verosimilitud parcial de Cox (Cox (1972), Martinussen y Scheike (2006)), 4.5. M´ ETODOS DE DIAGN ´ OSTICO EN EL MODELO DE COX 39 Los FP’s aproximan cada funci´on desconocida gpor una combinaci´on lineal de M polinomios xpj,j= 1,...,M. En el polinomio ordinario las potencias pjson restringidas para valores enteros positivos, mientras en modelamiento con FP’s se trabaja con valores positivos, nopositivos y valores fraccionales para pj. Un t´ıpico conjunto de potencias admisibles es dado por pj∈(−2,−1,−0,5,0,5,1,2,3), donde x0denota ln(x). M´as formalmente, un FP de grado Mes definido como FPM(x) = M X j=1 βjhj(x), donde β1,...,βMson los coeficientes de regresi´on y hjes recursivamente definido como h0(x) = 1, y hj(x) = (xpj,si pj=pj−1 hj−1(x)ln(x),si pj6=pj−1 Para M= 2 y pj6=pj−1:F P2(x) = β1xp1+β2xp2. Para M= 2 y pj=pj−1:F P2(x) = β1xp1+β2xp1ln(x). Programas para ajustar el modelo aditivo basado en FPs es evaluable en la plataforma de computaci´on estad´ıstica STATA (funci´on mfp), SAS (macro mfp8) y R (funci´on fp del paquete mfp). La implementaci´on en R se restringe para FP’s de grado 2 (Sauerbrei et al., 2006). 4.5. M´etodos de diagn´ostico en el modelo de Cox En un modelo de regresi´on lineal es f´acil definir un residuo. Sin embargo, en el modelo de regresi´on para datos de supervivencia la definici´on del residuo no es tan clara. Una serie de residuos se han propuesto para el modelo de Cox, que son ´utiles para examinar los diferentes aspectos del modelo (Klein y Moeschberger (1997), Therneau y Grambsch (2000)). En el modelo Cox (11), las restricciones naturales suponen verificar los siguientes aspectos: El logaritmo de la raz´on de riesgo no depende del tiempo, ln(λ(t, x)/λ0(t)) = β1X1+...+βpXp(riesgo constante). El riesgo de un individuo es proporcional al riesgo de otro individuo (riesgo proporcional). El logaritmo de la raz´on de riesgo y las covariables se relacionan linealmente (forma funcional lineal). 40 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES Una medida para evaluar la suposici´on de riesgos proporcionales puede ser realizada mediante m´etodos num´ericos o aproximaciones gr´aficas. Aqu´ı describimos brevemente los procedimientos basados en la gr´afica de los residuos de Schoenfeld escalado (Grambsch y Therneau, 1994), el estad´ıstico de prueba de no proporcionalidad de Therneau y Grambsch (2000) y los residuos martingala. 4.5.1. Verificaci´on del supuesto de riesgos proporcionales. 4.5.1.1. M´etodo gr´afico basado en los residuos de Schoenfeld escalado. Schoenfeld (1982) propone unos residuos para verificar la suposici´on de riesgo proporcional. Estos residuos son conocidos como residuos de Schoenfeld. El residuo de Schoenfeld es la diferencia entre el valor observado y esperado de la covariable en momento del tiempo (Therneau y Grambsch 2000). Sea el modelo de Cox extendido con efectos que var´ıan en el tiempo (Extended Cox model with time-varying coefficients), λ(t) = Y(t)λ0(t) exp(X0(t)β(t)) (28) donde β(t) es el efecto tiempo dependiente, que cuando no es constante, el impacto de una o m´as covariables en el riesgo puede variar sobre el tiempo. Pero la restricci´on β(t) = βimplica riesgo proporcional, por tanto, la gr´afica de β(t) versus el tiempo ser´ıa una l´ınea horizontal. Los residuos de Schoenfeld permiten detectar la variaci´on en el tiempo para un predictor de inter´es. En ausencia de empates estos residuos son iguales a la diferencia entre el vector de covariables observados y esperados para un evento en el tiempo tk(k= 1,...,d), erk=e Xk−E(e Xk/Rk), donde E(e Xk/Rk) = PlRke Xlexp(e βXl) PlRkexp(e βXl. En la presencia de pcovariables, los residuos de Schoenfeld erforman una matriz d×p, donde cada covariable ptiene un coeficiente estimado para cada evento del tiempo, βkp. Los residuos de Schoenfeld escalados (rs) se definen como el producto de la inversa del estimador de la matriz de varianza - covarianza del k-´esimo residuo de Schoenfeld y el k-´esimo residuo de Schoenfeld. Grambsch y Therneau (2000) muestran que E(rskp) + ˆ βp≈βp(tk), donde rskson los residuos de Schoenfeld escalados y ˆ βes un coeficiente estimado del modelo de Cox. Esto sugiere graficar rsk+ˆ βpversus tiempo, o alguna funci´on de tiempo g(t), como un m´etodo para visualizar la forma funcional de la variaci´on en el tiempo de los residuos Schoenfeld escalados para una covariable espec´ıfica. En el supuesto de 4.5. M´ ETODOS DE DIAGN ´ OSTICO EN EL MODELO DE COX 41 5e−02 5e−01 5e+00 5e+01 −0.10 −0.05 0.00 0.05 0.10 0.15 Time Beta(t) for edad 5e−02 5e−01 5e+00 5e+01 −1.5 −1.0 −0.5 0.0 0.5 Time Beta(t) for hb Gr´ afica 5. Residuos Schoenfeld escalado vs. g(t) = log(t), para la edad y Hb. riesgos proporcionales, los residuos se distribuyen alrededor de la l´ınea horizontal para un coeficiente βpconstante. Para facilitar la interpretaci´on de estos gr´aficos se superpone una curva de ajuste, utilizando alguna funci´on de ajuste local como lowess o loess (Cleveland, 1981). Si se cumple la hip´otesis de riesgo proporcional, los residuos deber´ıan agruparse de forma aleatoria a ambos lados del valor 0 del eje Y, y la curva ajustada deber´ıa ser pr´oxima a una l´ınea recta. En la gr´afica 5, se muestra la relaci´on entre rsk+ˆ βpyg(t) = log(t), para los datos de la edad y la hemoglobina. Las l´ıneas negras corresponden a la curva de ajuste de los residuos Schoenfeld escalado ±error est´andar aproximado mediante lowess y la l´ınea roja es la recta horizontal en el punto 0 del eje de las ordenadas. En esta gr´afica se observa que el efecto de la edad no var´ıa en el tiempo, ya que los residuos se aproximan a la l´ınea horizontal en el punto 0 del eje Y, lo que significa el cumplimiento de riesgo constante. En cambio el efecto de la hemoglobina disminuye y despu´es se incrementa a lo largo del tiempo, lo que contradice la suposici´on de riesgo constante a lo largo del tiempo de un modelo de Cox correctamente especificado. 4.5.1.2. Test de no-proporcionalidad de Therneau y Grambsch. Siguiendo la aproximaci´on gr´afica, Grambsh y Therneau (1994) introducen una versi´on del test del score basado en los residuos de Schoenfeld escalados. Escrito β(t) como una funci´on de regresi´on en g(t), los coeficientes tiempo dependientes del modelo de Cox extendido (28) pueden ser escritos como, βj(t) = βj+θj(gj−¯gj), j = 1,...,p (29) 42 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES donde ¯gjes la media de la gj(funci´on de tiempo especificado). Una t´ıpica aplicaci´on de este tipo de prueba es para gj=log(t). En la expresi´on (29) el inter´es es realizar una prueba de hip´otesis sobre la hip´otesis de riesgo proporcional global H0:θ= 0 y para una covariable espec´ıfica H0j:θj= 0, j= 1,...,p. Para realizar estas pruebas de hip´otesis, Therneau y Grambsch (2000) introducen dos estad´ısticos de prueba, una para la global y otro para una covariable especifica. Para la hip´otesis global, H0:θ= 0, el estad´ıstico de prueba de riesgo proporcional para todas las pcovariables es T=(g−¯g)0S∗IS∗0(g−¯g) dP(gk−¯g)2, donde S∗es la matriz de los residuos de Schoenfeld escalados, Ies la matriz de informaci´on y des eventos de tiempo. Para la hip´otesis de una covariable espec´ıfica H0j:θj= 0, el estad´ıstico de prueba de riesgo proporcional es Tj=(P(gk−¯g)rskj)2 dIjj P(gk−¯g)2, donde Ijj es el elemento de la matriz de informaci´on para la j-´esima covariable ydson los eventos de tiempo. Este estad´ıstico se distribuye asint´oticamente como una χ2 1. Por otro lado, sean los coeficiente de regresi´on tiempo-dependiente del modelo de Cox extendido (28) que puede ser escrito como βp(t) = βp+θpgp(t),(30) donde gp(t) es una funci´on de tiempo especificado previamente. Cuando las g0sson funciones conocidas, entonces el modelo con coeficientes (30) es a´un el modelo de Cox. Por tanto, el estimador de los par´ametros se puede obtener maximizando la funci´on de verosimilitud partial y las pruebas de hip´otesis se pueden realizar sobre las componentes tiempo-dependientes utilizando el test de score, la prueba de raz´on de verosimilitud o la prueba de Wald. (Martinussen y Scheike (2006), Therneau y Grambsch (2000)). Sea U= (U0 1, U0 2) la funci´on de puntaje, donde la primera componente es la derivada de la verosimilitud parcial con respecto a βy la segunda componente respecto a θ. Sea Ikl (k, l = 1,2) la matriz de informaci´on emp´ırica definida como una matriz de bloques reflejando dos vectores de par´ametros, . Para realizar la prueba de hip´otesis sobre H02 :θ= 0, θ= (θ1,...,θp) de los par´ametros de la expresi´on (30), se usa el estad´ıstico de prueba global (score test) para los efectos tiempo-dependientes definida como 4.5. M´ ETODOS DE DIAGN ´ OSTICO EN EL MODELO DE COX 43 T(G) = U0 2(ˆ βp,0)I−1 22 (ˆ β, 0)U2(ˆ βp,0), donde ˆ βdenota el estimador de m´axima verosimilitud parcial. Este estad´ıstico se distribuye como una χ2con pgrados de libertad bajo la hip´otesis nula. Para los datos de la gr´afica de los residuos de Schoenfeld escalados vs. logaritmo del tiempo (Gr´afica 5), los resultados (valores p) del estad´ıstico de prueba de no proporcionalidad de Therneau y Grambsch indican que el efecto de la edad en la funci´on de riesgo es constante (p = 0.359) y en cambio el efecto de la hemoglobina es de riesgo no constante (p = 0.009). En el package R, la funci´on cox.zph permite realizar la gr´afica de los residuos de Schoenfeld escalados versus una transformaci´on de la funci´on tiempo (g(t)) y obtener los estad´ısticos de prueba individual y global. Las transformaci´on de la funci´on de tiempo disponibles son la identidad g(t) = t,g(t) = log(t), rangos de eventos de tiempo y por defecto es 1-KM (KM: es el estimador de Kaplan-Meier). 4.5.2. Verificaci´on de la forma funcional lineal. En el modelo de Cox, uno de los residuos muy usuales para verificar la forma funcional del efecto de las covariables en la funci´on de riesgo son los residuos basados en martingalas. Barlow y Prentice (1988) provee el marco b´asico de los residuos martingala y el posterior trabajo de Therneau, Grambsch y Fleming (1990). Sea Ni(t) el proceso de conteo (n´umero de eventos observados) y la funci´on de intensidad acumulada Λi(t) = Rt 0λsds, con informaci´on adicional en t´erminos de p covariables Xi. Los residuos martingala Mi(t) se definen como la diferencia entre los procesos de conteo observados y esperados, Mi=Ni(t)−Ei(t), donde Ei(t) = Λi(t). Bajo la funci´on de intensidad de la forma (17) y usando los estimadores del modelo de Cox se pueden estimar los residuos martingala, Mi(t). Por otro lado, los residuos martingala se pueden construir bas´andose a su vez en los denominados residuos de Cox-Snell (rci), crmi=δi−brci, donde δies 1 si ocurre el evento, 0 caso contrario y los residuos de Cox-Snell se definen como el estimador de la funci´on de intensidad acumulado, brci=ˆ Λ(ti), i= 1,...,n. Si la muestra es grande, la suma de los residuos martingala es cero, son no correlacionados y el valor esperado es cero. Sin embargo, no se distribuyen de forma sim´etrica en torno a cero, aunque el modelo sea correcto, lo que dificulta la interpretaci´on de los gr´aficos. La gr´afica de los residuos martingala versus la covariable, bajo el supuesto de efecto lineal en el modelo, deben verificar que los residuos se distribuyen alrededor de un punto del eje y, sin que sugiera una curva de ajuste de forma funcional no lineal. En el gr´afica 6 se muestra los residuos martingala versus la edad y hemoglobina. Las l´ıneas negras corresponden a la curva de ajuste de los residuos aproximado mediante lowess y la l´ınea roja es la recta horizontal en el punto 0 del eje y. En esta gr´afica se observa que el efecto de la edad y de la hemoglobina no es lineal; ya que la curva de ajuste de los residuos versus la edad y hemoglobina (l´ıneas rojas) 44 4. MODELO DE REGRESI ´ ON DE COX CON M´ ETODOS FLEXIBLES 20 40 60 80 −1.0 −0.5 0.0 0.5 1.0 edad Residuo martingala 5 10 15 −1.0 −0.5 0.0 0.5 1.0 Hb Residuo martingala Gr´ afica 6. Residuos martingala vs. edad y Hb. no se aproximan a una l´ınea horizontal. Los cuales, significan que el efecto de las covariables en la funci´on de riesgo presentan una forma funcional no lineal. Por otro lado, Lin et al. (1993) y Wei (1984) sugieren una importante clase de estad´ısticos de prueba basado en la suma acumulada de los residuos martingala. Estos estad´ısticos son dise˜nados para investigar las diferentes salidas del modelo, incluyendo errores de especificaci´on de la funci´on de enlace y la forma funcional de las covariables (Martinussen y Scheike, 2006). Las martingalas bajo el supuesto de riesgo proporcional del modelo de Cox se pueden escribir como Mi(t) = Ni(t)−Zt 0 Yi(s) exp(X0 iβ)dΛ0(s). En el cual, usando los estimadores del modelo de Cox se puede estimar Mi(t) como ˆ Mi(t) = Ni(t)−Zt 0 Yi(s) exp(X0 iˆ β)dˆ Λ0(s). La idea ahora es mirar las diferentes funcionales de estos residuos estimados y ver si se comportan como deber´ıa bajo el modelo propuesto. Lin et al. (1993) define un proceso de residuales acumulativo bi-dimensional como Mc(t, z) = Zt 0 Kt z(s)dˆ M(s), donde Kz(t) es una matriz n×1 con elementos I(Xi1≤z) para i= 1,...,n, centr´andose aqu´ı en la primera covariable continua X1. En este caso, los residuos martingala son agrupados de forma acumulativa respecto al tiempo de seguimiento y valores de la covariable. Para resumir ´estas se puede integrar sobre el periodo de 4.5. M´ ETODOS DE DIAGN ´ OSTICO EN EL MODELO DE COX 45 tiempo y conseguir un proceso ´unicamente en z,Mc(z) = Rt 0K0 z(t)dˆ M(t), el cual puede ser graficado contra z. Para evaluar el proceso observado como inusual bajo el modelo propuesto, se puede graficar este a lo largo del tiempo como una realizaci´on bajo el modelo. Para mejorar a´un m´as la objetividad de la t´ecnica gr´afica, se puede completar con con un estad´ıstico de prueba llamada test de supremo de Mc(t); el cual, mide el extremo del proceso observado. Un valor demasiado grande de este test sugiere que la forma funcional lineal del efecto de la covariable es inapropiada, lo que significa que las variables no podr´ıan entrar en el modelo en la escala original; por tanto, esta variable requiere alg´un tipo de trasformaci´on para ser incluido en el modelo. En el package R, la funci´on cox.aalen permite realizar la simulaci´on de los residuos y obtener la gr´afica de los residuos acumulados vs. las covariables continuas, as´ı como el test supremo. Cap´ıtulo 5 Factores pron´osticos en LNH 5.1. Descripci´on de los datos En este trabajo se analizan los datos de 2160 pacientes mayores o iguales a 14 a˜nos de edad con diagn´ostico de linfoma no Hodgkin (LNH) que fueron diagn´osticados y tratados en el Instituto Nacional de Enfermedades Neopl´asicas (INEN), Lima-Per´u, entre 1990 y 2002. El tratamiento que hab´ıan recibido los pacientes seg´un la pr´actica cl´ınica habitual (seg´un el protocolo de tratamiento) fue generalmente quimioterapia en la mayor´ıa de los casos (91.2 %) y los restantes (8.8 %) hab´ıan recibido radioterapia y/o cirug´ıa. El esquema de quimioterapia fue generalmente (81.6 %) CHOP (ciclofosfamida, doxorubicina, vincristina y prednisona) y los restantes otros esquemas de quimioterapia. Las siguientes caracter´ısticas cl´ınicas (covariables) de los pacientes, documentados al diagn´ostico, fueron incluidos en el an´alisis: la edad (en a˜nos), g´enero (femenino, masculino), estado funcional (zubrod: 0, 1, 2, 3, 4), foco primario (primario: enfermedad ganglionar, extraganglionar), estadio cl´ınico (EC: I, II, III, IV), s´ıntomas B (fiebre, sudoraci´on nocturna o baja de peso, sin causa alguna), hemoglobina (Hb: en g/dl), leucocitos (leuco: en mil/mm3), linfocitos (linf: en %), y deshidrogenasa l´actica (DHL: en UI/L). No se incluye la β-2 microglobulina (β2M: en mg/L), debido a que la mayor´ıa de los pacientes no ten´ıan informaci´on al respecto. El tiempo de supervivencia (en meses) que es la variable a modelar en t´erminos de las covariables, fue calculado desde la fecha de diagn´ostico hasta la fecha de muerte o fecha del ´ultimo control que fue registrada en la historia cl´ınica. Los pacientes fallecidos se consideraron como eventos (no censurados) y los restantes como censurados. De los 2160 pacientes con LNH, 709 (32.8 %) pacientes hab´ıan fallecido y la mediana de seguimiento de los pacientes restantes fue de 12.6 meses. La mediana de supervivencia fue de 61.8 meses (IC95: 49.9 - 73.7) y la tasa de supervivencia a 5 y 10 a˜nos de 51.2 % y 41.7 % respectivamente. En la Gr´afica 7 se muestra la curva de supervivencia estimada mediante el m´etodo de Kaplan-Meier y la funci´on de riesgo acumulado. 47 48 5. FACTORES PRON ´ OSTICOS EN LNH 0 50 100 150 0.0 0.2 0.4 0.6 0.8 1.0 Curva de supervivencia global Meses Survival function 0 50 100 150 0.0 0.2 0.4 0.6 0.8 Función de riesgo acumulado Meses Cumulative Hazard Gr´ afica 7. Curva de supervivencia y riesgo acumulado de los pacientes con LNH. 5.2. Aplicando el modelo de Cox cl´asico En la Tabla 2 se muestran los resultados de la aplicaci´on del modelo de Cox, incluy´endose las siguientes variables: edad, g´enero (femenino, masculino), zubrod (0-1, 2-4), primario (ganglionar, extraganglionar), estadio cl´ınico (I-II, III-IV), s´ıntomas (A, B), Hb, ln(leucocitos), linfocitos y ln(DHL). Tabla 2. Resultados del modelo de Cox cl´asico. Variables ˆ βEE(ˆ β) Z p HR (IC95 %) Edad (a˜nos) 0.005 0.002 2.210 0.027 1.01 (1.00, 1.01) G´enero masculino 0.202 0.078 2.590 0.010 1.22 (1.05, 1.43) Zubrod 2-4 0.689 0.085 8.14 <0.001 1.99 (1.69, 2.35) Primario ganglionar -0.046 0.085 -0.550 0.580 0.96 (0.81, 1.13) Estadio cl´ınico III-IV 0.440 0.086 5.140 <0.001 1.55 (1.31, 1.84) S´ıntomas B 0.164 0.082 2.00 0.045 1.18 (1.00, 1.38) Hb -0.051 0.016 -3.16 0.002 0.95 (0.92, 0.98) ln(leucocitos) 0.211 0.069 3.03 0.002 1.24 (1.08, 1.42) Linfocitos -0.014 0.003 -4.32 <0.001 0.99 (0.98, 0.99) ln(DHL) 0.255 0.044 5.73 <0.001 1.29 (1.18, 1.41) LR test 351.00 AIC 9597.33 Nota: Las categor´ıas que no aparecen corresponden a las categor´ıas de referencia. Beta: Coeficiente de regresi´on. EE: error est´andar del coeficiente de regresi´on. HR: Hazard ratio. 5.3. APLICANDO EL MODELO DE COX CON P-SPLINES 55 Tabla 5. Prueba individual y global para la proporcionalidad, basado en el residual de Schoenfeld escalado del modelo de Cox con P-splines. La correlaci´on de Pearson es entre los residuos y g(t) para cada covariable. Variables correlaci´on de Pearson (rho) χ2valor-p pspline(edad, df = 4) 0.01664 0.43494 0.5096 G´enero masculino 0.00123 0.00220 0.9626 Zubrod 2-4 -0.05522 4.48729 0.0341 Primario ganglionar -0.01379 0.28662 0.5924 S´ıntomas B 0.02555 0.98132 0.3219 Estadio cl´ınico III-IV 0.01355 0.26061 0.6097 pspline(hb, df = 4) -0.00132 0.00263 0.9591 pspline(lleu, df = 4) 0.00780 0.10383 0.7473 pspline(linf, df = 4) -0.05655 4.63394 0.0313 pspline(ldhl, df = 4) -0.02377 0.72981 0.3929 Global NA 13.89722 0.1777 20 60 −4 −3 −2 −1 0 1 Edad Residuos martingala Fem Masc −4 −3 −2 −1 0 1 Género Residuos martingala 0−1 2−4 −4 −3 −2 −1 0 1 Escala ECOG Residuos martingala Ext Gang −4 −3 −2 −1 0 1 Primario Residuos martingala I−II III−IV −4 −3 −2 −1 0 1 EC Residuo martingala A B −4 −3 −2 −1 0 1 Síntomas Residuos martingala 5 10 15 −4 −3 −2 −1 0 1 Hb Residuos martingala 7 9 11 −4 −3 −2 −1 0 1 Ln(Leuco) Residuos martingala 0 40 80 −4 −3 −2 −1 0 1 Linfocitos Residuos martingala 3579 −4 −3 −2 −1 0 1 Ln(DHL) Residuos martingala Gr´ afica 13. Residuos martingala del modelo de Cox con P-spline. 56 5. FACTORES PRON ´ OSTICOS EN LNH 5.4. Aplicando el modelo de Cox con PF En la Tabla 6 se muestra los resultados del modelo de Cox con polinomio fraccional para aproximar el efecto de las covariables continuas (edad, Hb, ln(leuco), linfocitos y ln(DHL) en la funci´on de riesgo, incluy´endose las covariables categ´oricas como el g´enero, zubrod, primario, estadio cl´ınico y s´ıntomas. Tabla 6. Resultados del modelo de Cox con polinomio fraccional. Variables ˆ βEE(ˆ β) Z p HR (IC95 %) Edad: I((edad/100)3) 0.765 0.248 3.086 0.002 2.15 (1.32, 3.49) G´enero masculino 0.227 0.078 2.926 0.003 1.26 (1.08, 1.46) Zubrod 2-4 0.631 0.084 7.474 <0.001 1.88 (1.59, 2.22) Primario ganglionar -0.087 0.085 -1.023 0.306 0.92 (0.78, 1.08) EC III-IV 0.402 0.086 4.683 <0.001 1.50 (1.26, 1.77) S´ıntomas B 0.110 0.082 1.340 0.180 1.12 (0.95, 1.31) Hb1: I((Hb/10)−2) 0.675 0.172 3.915 <0.001 1.96 (1.40, 2.75) Hb2: I((Hb/10)−2×ln(Hb/10)) 0.649 0.164 3.955 <0.001 1.91 (1.39, 2.64) Ln(leuco): I((lleuco/10)1) 1.159 0.678 1.709 0.087 3.19 (0.84, 12.04) Linfocitos1: I((linf/10)0,5) -1.229 0.140 -8.758 <0.001 0.29 (0.22, 0.39) Linfocitos2: I((linf/10)2) 0.039 0.006 6.532 <0.001 1.04 (1.03, 1.05) Ln(DHL): I(ldhl/10)1) 2.600 0.451 5.762 <0.001 13.46 (5.56, 32.58) LR test 427.20 AIC 9525.13 Nota: Las categor´ıas que no aparecen corresponden a las categor´ıas de referencia. Beta: Coeficiente de regresi´on. EE: error est´andar del coeficiente de regresi´on. HR: Hazard ratio. En los resultados de la Tabla 6 se observa que los factores pron´osticos con efecto significativo (p <0.05) en la supervivencia global de los pacientes con LNH son casi todas las variables, a excepci´on del primario (p = 0.306), s´ıntomas (p = 0.180) y los leucocito (p=0.087). Los resultados de este modelo son similares al modelo con P-splines a excepci´on del efecto de los leucocitos. Para las covariables categ´oricas la tasa de riesgo de estas variables implican que, los pacientes de sexo masculino presentan un riesgo de mortalidad de HR=1.3 (IC5 %: 1.1-1.5) veces mas que los pacientes de sexo femenino. Los pacientes con zubrod 2-4 presentan un riesgo de mortalidad de HR=1.9 (IC95 %: 1.6-2.2) veces m´as que los pacientes con zubrod 0-1. Los pacientes con enfermedad avanzada (EC III-IV) presentan un riesgo de mortalidad de HR=1.5 (IC95 %: 1.3-1.8) veces m´as que los pacientes con enfermedad temprana (EC I-II). Seg´un el m´etodo de polinomio fraccional cada covariable continua requiri´o un tipo de trasformaci´on. La edad fue dividida por 10 y despu´es requiri´o una trasformaci´on de potencia c´ubica, la Hb fue transformada en (I(Hb/10))−2+I((Hb/10))−2× ln(Hb/10), el logaritmo del n´umero de leucocitos fue dividida por 10, el porcentaje de linfocitos fue transformada en I((linfocitos/10))0,5+I((linfocitos/10))2y el logaritmo de DHL fue dividida por 10. 5.4. APLICANDO EL MODELO DE COX CON PF 57 En la Gr´afica 14 se muestra la forma funcional de las covariables en la raz´on de riesgo. Para la edad, el riesgo de mortalidad en menores de 60 a˜nos es menor a HR=1, sin embargo, el riesgo se incrementa ligeramente despu´es de los 60 a˜nos de edad. Para la Hb ≥12g/dl, el riesgo de mortalidad es aproximadamente menor a HR=1, sin embargo, el riesgo se incrementa cuando la Hb disminuye despu´es de 12g/dl. Para el n´umero de leucocitos, el riesgo de mortalidad es aproximadamente menor a HR=1 para menores a 10 mil, y cuando el n´umero de leucocitos se incrementa despu´es de los 10 mil leucocitos el riesgo de mortalidad tambi´en se incrementa (efecto de hiperleucocitosis); aunque aqu´ı no se detecta el riesgo de mortalidad por leucopenia como fue identificado por los splines. Para el porcentaje de linfocitos, el riesgo de mortalidad se incrementa cuando el porcentaje de leucocitos disminuyen despu´es 20 %, as´ı mismo, se incrementan despu´es de los 60 %. Para la deshidrogenasa l´actica (DHL), el riesgo de mortalidad para un DHL menor a 240UI/L es menor a HR=1, sin embargo, cuando el nivel de DHL se incrementa despu´es de 240UI/L el riesgo tambi´en se incrementa. 20 60 0 1 2 3 4 5 Edad Hazar Ratio 1.0 1.4 1.8 0 1 2 3 4 5 género Hazar Ratio 1.0 1.4 1.8 0 1 2 3 4 5 zubrod Hazar Ratio 1.0 1.4 1.8 0 1 2 3 4 5 primario Hazar Ratio 1.0 1.4 1.8 0 1 2 3 4 5 EC Hazar Ratio 1.0 1.4 1.8 0 1 2 3 4 5 síntomas Hazar Ratio 5 10 15 0 1 2 3 4 5 Hb Hazar Ratio 7 9 11 0 1 2 3 4 5 ln(leuco) Hazard Ratio 0 40 80 0 1 2 3 4 5 Linfocitos Hazard Ratio 3 5 7 9 0 1 2 3 4 5 ln(DHL) Hazar Ratio Gr´ afica 14. Forma funcional de la raz´on de riesgo mediante el modelo de Cox con PF. En el Gr´afico 15 se muestra los residuos martingala para el modelo de Cox con polinomio fraccional. Los resultados indican que los residuos para cada covariable son aproximadamente constantes, lo cual verifica que los efectos de las covariables son adecuadamente aproximados utilizando el m´etodo de polinomio fraccional. 58 5. FACTORES PRON ´ OSTICOS EN LNH 20 60 −5 −4 −3 −2 −1 0 1 Edad Residuos martingala Fem Masc −5 −4 −3 −2 −1 0 1 Género Residuos martingala 0−1 2−4 −5 −4 −3 −2 −1 0 1 Escala ECOG Residuos martingala Ext Gang −5 −4 −3 −2 −1 0 1 Primario Residuos martingala I−II III−IV −5 −4 −3 −2 −1 0 1 EC Residuos martingala A B −5 −4 −3 −2 −1 0 1 Síntomas Residuos martingala 5 10 15 −5 −4 −3 −2 −1 0 1 Hb Residuos martingala 7 9 11 −5 −4 −3 −2 −1 0 1 Ln(Leuco) Residuos martingala 0 40 80 −5 −4 −3 −2 −1 0 1 Linfocitos Residuos martingala 3579 −5 −4 −3 −2 −1 0 1 Ln(DHL) Residuos martingala Gr´ afica 15. Residuos martingala del modelo de Cox con PF. 5.5. Comparaci´on de los modelos En la Tabla 7 se muestran los resultados resumidos del ajuste de los tres modelos (Cox cl´asico, Cox con P-splines y Cox con polinomio fraccional). Aqu´ı se comparan los factores pron´osticos identificados, la forma funcional del efecto de las covariables en la funci´on de riesgo y la selecci´on del mejor modelo mediante AIC (criterio de informaci´on de Akaike) Los factores pron´osticos identificados mediante el modelo de Cox cl´asico fueron casi todas las covariables incluidas en el an´alisis, a excepci´on del foco primario (p=0.583); sin embargo, en este modelo las variables como zubrod, primario, s´ıntomas, ln(leucocitos) y ln(DHL) no cumpl´ıan los supuestos de riesgo proporcional. As´ı mismo, los residuos martingala muestran que el efecto de estas variables en la supervivencia presentan una forma funcional no lineal, lo que suger´ıa el uso de m´etodos m´as flexibles para aproximar el efecto de las covariables en la funci´on de riesgo. En cambio con el modelo de Cox con P-splines y modelo de Cox con polinomio fraccional las covariables con efecto significativo fueron casi todas a excepci´on del primario y s´ıntomas (p >0.05). Si bien ambos modelos describen la forma funcional 5.5. COMPARACI ´ ON DE LOS MODELOS 59 no lineal de los efectos de las covariables, el ajuste con P-splines describe mejor la forma funcional en algunas covariables que no es identificado por el polinomio fraccional. Un ejemplo, de esto es para logaritmo de leucocitos donde P-splines describe dos grupos de riesgo (para <3mil y >10mil) y el polinomio fraccional solo identifica un grupo de riesgo linealmente creciente para leucocitos mayores de 10mil. Los residuos martingala para ambos modelos (modelo de Cox con P-splines y modelo de Cox con polinomio fraccional) seg´un las covariables no muestran un patr´on que sugiera la existencia de alguna forma funcional que evidencie la falta de ajuste de los datos. Para ambos modelos los residuos martingala por cada covariable son aproximadamente constantes en consecuencia los ajustes de ambos modelos son adecuados. Sin embargo, de acuerdo al criterio de informaci´on de Akaike (AIC), el modelo de Cox con splines presenta ligeramente una menor AIC (9523.97) en comparaci´on al modelo de Cox con polinomio fraccional (AIC: 9525.13). Finalmente en la Tabla (8) se muestra los factores pron´osticos para la supervivencia de los pacientes con Linfoma no Hadgkin, bajo el modelo de Cox con P-splines. Las variables con efecto significativo para la supervivencia a un nivel de significaci´on de significaci´on de 5 % (α= 0,05) fueron: la edad (p = 0.028), g´enero (p = 0.009), zubrod (p <0.001), EC (p <0.001), Hb (p = 0.006), leucocitos (p = 0.038) y linfocitos (p <0.001), as´ı mismo, la DHL (p = 0.089) por ser cl´ınicamente relevante por su significado pron´ostico en los LNH. 60 5. FACTORES PRON ´ OSTICOS EN LNH Tabla 7. Resultados del modelo de Cox con polinomio fraccional. Modelo de Cox Modelo de Cox con P-spline Modelo de Cox con PF Variables sin transformaci´on pHR P-spline pHR PF pHR Edad: √0.027 1.01 lineal 0.007 I((edad/100)3) 0.002 2.15 - - - no-lineal 0.028 pspline( , df=4) - - - G´enero masculino √0.009 1.22 √0.009 1.23 √0.003 1.26 Zubrod 2-4 √<0.001 1.99 √<0.001 1.86 √<0.001 1.88 Primario ganglionar √0.583 0.95 √0.109 0.89 √0.306 0.92 EC III-IV √<0.001 1.55 √<0.001 1.49 √<0.001 1.50 S´ıntomas B √<0.001 1.18 √0.140 1.13 √0.180 1.12 Hb: √0.002 0.95 lineal 0.082 I((Hb/10)2)<0.001 1.96 - - - no-lineal 0.006 pspline( , df=4) I((Hb/10)−2×ln(Hb/10)) <0.001 1.91 ln(Leucocitos): √0.002 1.24 lineal 0.082 I(ln(leuco)/10) 0.087 3.19 - - - no-lineal 0.038 pspline( , df=4) - - - Linfocitos: √<0.001 0.99 lineal <0.001 I((linf/10)0,5)<0.001 0.29 - - - no-lineal <0.001 pspline( , df=4) I((linf/10)2)<0.001 1.04 ln(DHL): √<0.001 1.29 lineal <0.001 I(ln(DHL/10)) <0.001 13.46 - - - no-lineal 0.089 pspline( , df=4) - - - LR test 351.00 455.00 427.2 AIC 9597.33 9523.97 9525.13 Nota: Las categor´ıas: ln:logaritmo natural,√: sin trasformac´ı´on. 5.5. COMPARACI ´ ON DE LOS MODELOS 61 Tabla 8. Factores pron´osticos para LNH, bajo el modelo de Cox con P-spline. Variables p HR (IC95 %) Edad 0.028 P-splines ( , df=4) G´enero masculino 0.009 1.23 (1.05, 1.44) Zubrod 2-4 <0.001 1.86 (1.57, 2.20) Primario ganglionar 0.190 0.89 (0.75, 1.06) EC III-IV <0.001 1.49 (1.26, 1.76) S´ıntomas B 0.140 1.13 (0.96, 1.33) Hb 0.006 p-splines ( , df=4) ln(leuco) 0.038 p-splines ( , df=4) Linfocitos <0.001 p-splines ( , df=4) ln(DHL) 0.089 p-splines ( , df=4) LR test 455.00 AIC 9523.97 Nota: Las categor´ıas que no aparecen corresponden a las categor´ıas de referencia. HR: Hazard ratio. IC95 %: Intervalo de confianza al 95 %. Cap´ıtulo 6 Discusi´on y conclusiones El amplio uso de los modelos tradicionales para el an´alisis de supervivencia ha contribuido al desarrollo de m´etodos m´as sofisticados desde t´ecnicas simples a m´as complejas, las cuales han crecido r´apidamente durante los ´ultimos a˜nos para un mejor modelamiento, facilitado por el r´apido desarrollo de la tecnolog´ıa computacional Si bien el modelo de Cox (1972) es una herramienta muy utilizada para determinar el efecto de las covariables en muchos contextos estad´ısticos, este modelo est´a sujeto al cumplimiento de los supuestos como son: riesgo proporcional, covariables invariantes en el tiempo y que la estructura de la relaci´on entre la funci´on de riesgo y las covariables sea lineal. Sin embargo, estas condiciones o restricciones no necesariamente se cumplen en muchas aplicaciones. En este sentido, la no-linealidad puede ser tan frecuente como el no cumplimiento de riesgos proporcionales; como algunos autores refieren uno puede ser consecuencia del otro, es decir, si no hay proporcionalidad es muy posible que tampoco haya linealidad (Keele, 2010). En consecuencia, si el supuesto de riesgos proporcionales no se cumple, el modelo de Cox cl´asico no es el m´as adecuado, entonces el modelo de Cox estratificado, modelo de Cox extendido con variable tiempo-dependiente, modelo de odds proporcional y modelo log-log´ıstico o modelo de Cox ponderado podr´ıan ser una alternativa, pero se debe tener en cuenta en todos estos modelos la forma lineal del efecto de las covariables. En este trabajo, se utilizaron m´etodos m´as flexibles como son: m´etodo de suavizamiento P-spline y polinomio fraccional, debido a que en nuestros datos, la forma funcional no satisface el supuesto de relaci´on lineal en el modelo de Cox cl´asico. En cambio utilizando el modelo de Cox con P-splines y el modelo de Cox con polinomio fraccional se obtuvieron una mejor aproximaci´on de los efectos de las covariables en la funci´on de riesgo. En consecuencia, la raz´on de riesgo para cada covariable continua presenta una estructura de relaci´on cuya forma funcional es no lineal. Los factores pron´osticos con efecto significativo para la supervivencia en LNH fueron: la edad, g´enero, zubrod, estadio cl´ınico (EC), nivel de hemoglobina (Hb), leucocitos, linfocitos y la deshidrogenasa l´actica (DHL) como en el modelo de Cox con P-splines y el modelo de Cox con polinomio fraccional. Los cuales, concuerdan con 63 64 6. DISCUSI ´ ON Y CONCLUSIONES los reportados en la literatura para esta patolog´ıa (Nicolaides, Dimos y Pavlidis, 1998; Rebasa, 2001) El modelo de Cox con P-splines y el modelo de Cox con polinomio fraccional aproximan bien el efecto de las covariables, sin embargo, el m´etodo basado en los P-splines describe mejor la forma funcional de la raz´on de riesgo para las covariables continuas que el modelo de Cox con polinomio fraccional, siendo esta una muy buena alternativa en situaciones donde la forma funcional es no lineal. Cabe resaltar que los puntos de corte (HR=1) determinados para las covariables continuas mediante estos m´etodos se aproximan a los puntos de corte definidos cl´ınicamente como grupos de peor pron´ostico. Seg´un los resultados del modelo de Cox con P-splines, los pacientes mayores de los 60 a˜nos de edad tienen un peor pron´ostico, el cual coincide con el punto de corte definido para clasificar a los pacientes seg´un la edad en grupos de mayor riesgo. Para la Hb baja (<12g/dl) y los valores elevados de la DHL (>240U/L) los puntos de corte obtenidos coinciden con los puntos de corte definidos cl´ınicamente para un peor pron´ostico. Sin embargo, para los valores de los leucocitos y los linfocitos existen dos puntos de corte que muestran un mayor riesgo de mortalidad: i) leucocitos menores de 3mil y mayores de 10mil, y ii) linfocitos menores de 20 % y >60 %. Estos grupos de pronostico deber´ıan ser considerados en la pr´actica cl´ınica al momento de clasificar a los pacientes en grupos de pron´ostico. Finalmente, este trabajo tiene algunas limitaciones en cuanto a la base datos disponible para realizar el an´alisis. Todos los datos fueron recopilados retrospectivamente de las historias cl´ınicas de los pacientes; en las cuales la mayor´ıa de los datos no fueron registrados de acuerdo a los objetivos de este estudio. Resultado de esto son los diferentes criterios de clasifcaci´on histopatol´ogica que no han permitido incluir en el an´alisis las variables como tipo histol´ogico, grados de agresividad e inmunofenotivo (tipo celular). As´ı mismo, datos como las β-2 microglobulinas que junto con tipo celular son factores pron´osticos en este grupo de pacientes. Como trabajo posterior desde el aspecto cl´ınico, se podr´ıa plantear realizar el an´alisis de los factores pron´osticos en los linfomas agresivos, principalmente linfomas de c´elulas grandes B difuso, que seg´un la clasifcaci´on de la Organizaci´on Mundial de la Salud (clasificaci´on actual de los LNH) representa el 80 % de los LNH. Desde el aspecto metodol´ogico se podr´ıa plantear realizar un an´alisis de los factores pron´osticos utilizando el modelo de Cox con splines penalizados para aproximar la forma funcional del efecto de las covariables, as´ı como aproximar los coeficientes tiempo-dependiente para las variables que no son constantes en el tiempo.