Usando R para analizar la distribución de los número primos

La Teoría de Números es sin duda, la reina de las matemáticas, aún a día de hoy existen numerosos enunciados muy fáciles de entender, pero a los que aún no hemos sido capaces de ofrecer una demostración definitiva.

En esta entrada quería hacer la presentación de una librería que existe en R y que se llama primes, que permite dar una visión un poco más geométrica a algunos problemas de esta maravillosa teoría, que debería ofrecerse tanto a nivel de facultad o inferior.

En este caso voy a hacer referencia al libro A Classical Introduction to Modern Number Theory, donde en sus primeras páginas se al conocido como Teorema de los Números Primos, demostrado al final del XIX por Hadamard y Vallé Poussin de modo independiente (nótese que estas cosas sucedían en esta época porque no existía internet …)

El teorema viene a establecer lo siguiente:

π(x) y x / ln(x) son infinitésimos equivalentes para x suficientemente grande

Donde π(x) es el número de primos entre 1 y x

Estoy seguro que nadie se sienta una tarde de verano y se pone a pensar sobre ¿Qué teorema demuestro hoy? Para llegar a un resultado como el anterior, seguro que pintaron en un plano cartesiano, por un lado los datos de la función π(x) y por otro los datos de la función x / ln(x) y se dieron cuenta que cuando x era muy grande, ambas funciones parecían que convergían a valores parecidos, pero claro, ponte tú a contar ahora, por ejemplo el número de primos entre por ejemplo 1 y 15.000

Para la desgracia de Hadamard y Vallé Poussin, el no tener en la época ordenadores y lenguajes como R que implementen funciones como las anteriores, daba además una dificultad extra al problema. Supongamos que se desea dibujar las gráficas de las funciones π(x) (en azul) y de f(x) = x / ln(x) (en naranja) en el intervalo 1 a 10.000.000 muestreando de 10.000 en 10.000, se podría hacer a partir de una función que existe en la librería que se denomina is_prime() generaría, con unas pocas líneas de código R una gráfica como la de la cabecera. Esta función indica, dada una lista de número cuál es primo y cuál no lo es y es mucho más recomendable su uso que otra función que dispone esta librería y que se denomina prime_count(), que en teoría estima el número de primos hasta una cierta cantidad, pero hace uso de los resultados existentes hasta día de hoy como el que se comenta en el siguiente párrafo, por lo que para la prueba que aquí se plantea no estaría demostrando nada.

Se observa que parece que las gráficas se separan, pero eso no ocurre indefinidamente, de hecho, se consigue acotar la función π(x), con una prueba de Pierre Durart en 2010 que consigue afinar un poco más estos resultados del siguiente modo:

x/(ln x – 1) < π(x) Si x > 5392

π(x) < x / (ln x – 1.1) Si x > 60183

A continuación se representa “nuestra” función π(x) , con las cotas anteriormente mencionadas bajo el mismo intervalo y muestreo de observaciones hasta obtener el siguiente gráfico:

Nótese que las 2 cotas, superior e inferior, están super-ajustadas de hecho, se hace una ampliación en la franja de los 8M y 9M se observa lo siguiente:

Un control prácticamente absoluto de la función π(x) que está restringida a las 2 franjas, mucho más exactas que el primer resultado de Hadamard y Vallé Poussin y es que la matemática moderna actual sigue aún obteniendo grandes resultados, donde antes no se podía.

Finalmente cabe comentar que esta función, de “contar” números primos y las desigualdades que se han ido demostrando, para poder conocer cómo se distribuyen estos, es otro de los grandes problemas de las matemática que lleva en boga los últimos 200 años y que sin duda, seguirá dando resultados aparentemente sencillos, pero para los que se requiere grandes dosis de ingenio.

Dado que existen este tipo de recursos en R (y creo que en python también), prometo comentar en próximos post, más resultados de esta índole que espero que os gusten.

Miradas al cielo de Madrid

Dado que hasta ahora Verano, ha habido pocos días sin nubes en la alegre Comunidad de Madrid, me he decidido, cámara en mano a tratar de hacer fotografía desde la misma ciudad de Madrid a pesar de su enorme contaminación lumínica. Un cielo que prácticamente no deja entrever estrellas de magnitud más allá de la 3. A pesar de todo insistí y esto es lo que obtuve:

La foto está tomada desde el suelo de mi terraza en Manoteras, donde tengo muchos edificios en frente, las correspondientes farolas de la zona para que la gente no tropiece y eso, y vecinos alrededor que enciende y apagan continuamente sus lucecitas, ...

Pues bien, si os fijáis en la foto veréis un punto verde casi en el centro, de hecho, si ampliais la foto veréis como ese punto verde crece y deja entrever como un anillo verde, ... Se trata de una nebulosa denominada M57 de magnitud 8.8 es decir, que ni de coña la podéis ver con vuestros ojitos.

En fotos buenas con más aumentos y no en Madrid, veríais tonalidades más rojizas y cálidas que las verdes, pero vamos, qué ibais a esperar de un cielo Bortle 7 como el que tenemos aquí.

Esta nebulosa está en la constelación de Lyra que es una constelación típica del verano cuya estrella más brillante Vega, es básicamente mi estrella preferida. Una constelación pequeña pero llena de objetos interesantes como este muy fácil de encontrar, de hecho, el punto verde este de m57 se encuentra en el extremo del pequeño trapecio que constituye la constelación y entre las estrellas beta y gamma de dicha constelación, muy fácil de encontrar con la cámara.

Para sacar esta foto me he valido de los mismos elementos que el anterior post de Astrofotografía: cámara reflex, montura star tracker, teleobjetivo de 300 mm sigma apo, y trípode, lo que sí he usado para limpiar de "caquita" la foto ha sido el programa pixinsight además se hacer un apilado sencillo de sólo 7 imágenes de 20 segundos de exposición, 5 darks, 10 flats y 10 bias. Con eso obtengo el resultado anterior, frente al resultado de una única foto con un mínimo tratamiento, que da una imagen triste tal y como se muestra a continuación:

Esta última visión capta la auténtica contaminación lumínica que no nos deja ver nada, una única toma e 20 segundos que si no se trata, ná de ná

Lo que no te contaron sobre la interpretabilidad de los Support Vector Machine

Uno de los modelos que más me fascinan, y no por su potencia predictiva, son los conocidos en español como Máquinas de Vectores Soportes. Lo me fascina de estos modelos son las ideas geométricas que subyacen en su entrenamiento y en el algoritmo que lo soportan, pero no voy a escribir sobre esto en esta entrada, ya que me quiero referir más a los outputs, es decir, qué handicaps ofrece la interpretación de un modelo de este tipo cuando se desea aplicar en problemas de datos.

Antes de nada hay que comentar que los SVM estuvieron muy de moda en kaggle hace más de un lustro porque ofrecían de un modo rápido predicciones muy potentes, de hecho, su origen se remonta a los años 90, en los que fueron usado por Vapnik para reconocimiento de escritura manual, ofreciendo un gran ajuste y ahorro de recursos frente a las Redes Neuronales que en aquella época estaban un poco de capa caída.

Cortes, C.; Vapnik, V. Support-vector networks Machine Learning, 20(3), 273-297

Mi primer acercamiento a estos modelos venían del "antiguo mundo de la minería de datos", donde ya por el año 2008, podían utilizarse en herramientas de minería de datos como SAS Enterprise Miner. En aquella época se abogaba por el uso de estos modelos en Riesgo de Crédito en aras de su potencia como alternativa los clásicos Credit Scoring, pero no llegó a prosperar la idea en el mundo bancario, donde los reguladores, en materia de riesgo, van con mucha prudencia a la hora de elegir modelos para temas donde se juega con enormes cantidades de dinero y no les falta razón, porque sin duda ellos si saben de las consecuencias de aplicar un mal modelo cuando se conceden préstamos, o que nadie pueda dar una explicación sólida cuando ocurran sucesos extraños derivados de los modelos que se aplican.

Bajo mi punto de vista, este tipo de modelos, en su versión más sencilla, ofrecen una falsa visión de interpretabilidad, que claro, cuando hay 2 variables como se muestra en la figura anterior, al final lo que se genera es, en los SVM de tipo lineal, una inecuación que separa las regiones del tipo a1 · x1 + a2 · x2 > 0, de modo que se puede saber qué, aparte de generarse una ecuación de ligadura entre todas las variables explicativas, el signo de los coeficientes a1 y a2 indica de si la variable favorece hacia una región u otra. Incluso en el caso de tener un número n de variables y crear fronteras del tipo: a1 · x1 + a2 · x2 + ... + an · xn> 0; se llegan a interpretaciones más o menos entendibles donde se tiene información no sólo del signo de cada una de las variables, sino incluso de grupos de variables entre sí, eso sí, sin darnos cuenta, hemos casi perdido el tema de la importancia de la variable aunque con un "poco de maña", se puede recuperar el concepto.

Pues bien, entonces ¿Qué problema ofrecen estos modelos? Estos modelos adolecen de lo que se conoce como Curse of Dimensionality de lo que se tratará en otro momento, pero ofrecen en mi opinión otro tipo de debilidades y es que en su versión lineal, no se va a poder hacer mucho para los casos generales de poblaciones no separables, es decir, en general los datos no van a estar agrupados en 2 subpoblaciones como en la figura anterior y se tendrá un lío enorme difícil de graficar. Incluso en situaciones de 2 variables, puede ocurrir que nos encontremos ante la necesidad de separar una población simulada como la siguiente:

En este caso, la versión lineal no dará ningún resultado satisfactorio, sin embargo es aquí donde los SVM a través del concepto de kernel y de lo que en geometría se conoce como "embedimiento", son capaces de generar fronteras capaces de curvarse cuando se proyectan en el plano y lograr así la separación, pero claro, el problema ahora está en ¿A qué precio se logra eso? El precio que hay que pagar es que no se va a obtener una ligadura claramente interpretable como en el caso anterior. Es decir, en el caso más general, donde alguien use un kernel genérico ante un número relativamente alto de variables 5 o más, prácticamente, no va a tener control sobre la frontera que crean los SVM y cómo dicha frontera podrá ser finalmente representada en el espacio original de 5 o más variables.

Ante una población como la de la anterior figura, en R se podría resolver del siguiente modo, mediante una búsqueda de un kernel radial (no lineal):

tune.out = tune(svm,factor(y)~.,data=xtrain,kernel="radial",sacale=FALSE,type='C-classification', ranges=list(cost=c(0.001, 0.01, 0.1, 1)))

Que lugar a la siguiente frontera, que gracias a que se está en dimensión 2, hasta podría dibujarse:

Sin embargo, se observa que ahora no hay 2 sino 3 regiones, con lo que se necesitan 2 inecuaciones posiblemente de tipo no lineal (lo cual complica aún más el tema), la interpretación además se complica incluso si se quiere plasmar algebraicamente la representación anterior. Se puede también buscar entre los vectores soporte que generan la frontera a ver si ofrecen algo de luz, pero en casos generales, no tendremos más que a lo sumo (y con suerte dependiendo de los datos) alguna que otra representación bonita de difícil interpretación.

En este sentido los SVM apelan a la continua lucha simplicidad vs interpretabilidad y seguramente en kaggle, los SVM de carácter no lineal habrían ganado concursos frente a regresiones logísticas sencillas, pero el problema está en ¿Quién entiende el modelo? ¿Hasta dónde hay que sacrificar la interpretabilidad en aras a la predictividad? Son preguntas sobre las que hay que tomar decisiones conjuntas entre los DS que aplican conocimientos y aterrizan conceptos complejos y la gente de negocio que son los que tienen la pasta y toman la decisión de aplicarlos o no.