Aunque el concepto de minería de datos esté casi indisolublemente asociado al de bases de datos enormes, en la práctica, el análisis y desarrollo de los modelos se realizan sobre muestras pequeñas.

Esencialmente, para lo que nos ocupa, es pequeño un conjunto de datos que cabe en la RAM de un PC. Actualmente son habituales las máquinas con 1 GB. A modo de comparación, la base de datos de clientes de una de las mayores compañías españolas y en la que trabajé hace un tiempo venía a ocupar 5 GB.

De acuerdo con la ley de Moore, dentro de [18 * log( 5 / 1 ) =] 41,8 meses ya podrá considerarse pequeña.

Pero durante 41,8 meses todavía tendremos que seguir muestreando datos.
Por eso, en este artículo se discuten procedimientos para realizar extracciones aleatorias mediante el muestreo aleatorio simple –la más básica de tales técnicas– de algunos de los gestores de bases de datos con los que es más probable encontrarse:

• DB2:

  select * from t tablesample bernuilly (p);

• Mysql:

  select * from t order by rand() limit n;

• Oracle:

  select * from t sample(p);

• PostgreSQL:

  select * from t order by random() limit n;

• SQL Server:

  select top n * from t order by newid();

En las queries anteriores, p representa un porcentaje –un valor entre 0 y 100– y n, un determinado número de filas.
El interesado en realizar muestreos de dicha manera, deberían considerar la información anterior como punto de partida para investigar, según proceda:

• las opciones adicionales de las que disponen DB2 u Oracle
• la validez estadística de los muestreos realizados con newid(), rand() o random()
• tener en cuenta que, a partir de la versión 2005, SQL Server soporta el comando tablesample.

De todos modos, los que no podemos salir de casa sin Python en el pendrive siempre tenemos la opción de volcar las tablas a un fichero de texto y recurrir a una versión del siguiente script:

import random
from sys import argv, exit
import os
if len(argv) == 1:
  print("Argumentos: path fichero umbral")
  exit(1)
random.seed()
path = argv[1]
nom_f_entrada = argv[2]
umbral = float(argv[3])
nom_f_salida = "muest_" + nom_f_entrada
f_entrada = open(os.path.join(path, nom_f_entrada), "r")
f_salida  = open(os.path.join(path, nom_f_salida), "w")
f_salida.write(f_entrada.readline())
l = f_entrada.readline()
while l:
	if (random.random() < umbral):
	f_salida.write(l)
	l = f_entrada.readline()
	f_entrada.close()
f_salida.close()

Nota: esta entrada forma parte de las que aparecían en un antiguo blog mío ya extinto cuyo contenido trato ahora de recuperar. Algo de cuanto en ella se lee, por lo tanto, puede que huela a rancio (como la que se refiere al tamaño de la memoria habitual en los ordenadores). Pero lo más debería todavía poder sostenerse en pie varios años después.

Tengo acceso a una máquina que, aunque anda un poco corta de memoria, cuenta con ocho CPUs. Tenía unas simulaciones bastante pesadas que correr y quise aprovechar su naturaleza perfectamente paralelizable. Y, de paso, hacer con R lo mismo por lo que he visto a un consultor de SAS cobrar a razón de 3.000 dólares diarios.

En el fondo, es una trivialidad. Supongamos que la función que implementa la simulación se llama foo. Habitualmente, haríamos

n.iter <- 1000
resultados <- replicate( n.iter, foo() )

y las simulaciones se ejecutarían secuencialmente. Con el paquete multicore, la tarea podría distribuirse por 6 (por ejemplo: es que si uso las ocho me riñen) de las CPUs así:

library( multicore )
resultados <- mclapply( 1:n.iter, foo, mc.set.seed = TRUE, mc.cores = 6 )

Notas:

• La función foo tiene que admitir un parámetro para poder usar mclapply: no existe una versión paralelizable de replicate en el paquete. Véase una posible solución a este problema en el ejemplo del final.
• La opción mc.cores limita el número de CPUs que se usarán en la paralelización.
• La opción mc.set.seed es muy importante: sin ella, cada subproceso compartirá la semilla y los valores aleatorios que generará serán los mismos.  ¡Normalmente no es eso lo que se quiere!

El paquete multicore incluye algunas funciones de más bajo nivel (parallel, fork, etc.) que podrían utilizarse para tareas más específicas. Pero mclapply es cómoda, simple y funciona estupendamente. Hay que hacer notar que el paquete está basado en la función fork tal cual está implementada en los sistemas operativos POSIX. Windows, lástima, no lo es.

Y para acabar, números:

library( multicore )
n.iter <- 1000
foo <- function(i) mean( rnorm( 100000 ) )  

system.time( resultados <- replicate( n.iter, foo() ) )
# 16.47 segundos en mi sistema (elapsed)

system.time( resultados <-
  mclapply( 1:n.iter, foo, mc.set.seed = TRUE, mc.cores = 6 ) )
# 4.34 segundos

RapidMiner es, posiblemente, la plataforma de minería de datos libre que mejor reputación goza. Hasta la publicación de la versión 5 le veía un pequeño problema: tenía una interfaz bastante poco intuitiva.

Hasta hace pocos días le veía otro: no podía extenderse —al menos de una manera obvia— programando en Java o, preferiblemente, R. Sin embargo, el módulo de integración de R con Rapidminer ya está listo y su lanzamiento va a ser el plato fuerte de RCOMM 2010, la conferencia de usuarios de Rapidminer (oficialmente, RapidMiner Community Meeting And Conference).

Estaremos muy atentos a ver qué da de sí esta integración tan esperada (de algunos, claro).

Si yo fuera responsable de una empresa de consultoría de, digamos (por mencionar un número redondo), 7.000 empleados, propondría la siguiente iniciativa:

• Crear un blog público de asuntos relacionados con la actividad de la empresa.
• Invitar a los empleados a mandar posibles entradas (originales, etc.).
• Cada día, publicar la mejor de las entradas recibidas.
• Recompensar a su autor con alguna chuche: un día de vacaciones, 100 euros, etc.

¿Ventajas? Supongo que muchas. De hecho, IBM patrocina una iniciativa similar en la que se publican artículos muy interesantes.

¿Inconvenientes? Todos, imagino, para los de las sinapsis cavernícolas.

¡A saber qué opinará Raúl de esto!

Una función muy útil de R es ifelse:

val <- 0
var <- ifelse( val == 1, "uno", "cero" )
print( var )

Un programador en SAS haría algo así como

%macro test(val);
 %if &val=1 %then %let var=one;
 %else %let var=zero;
 %put &var;
%mend;

%test(0);

SAS, sin embargo, recomienda hacerlo así:

%let val=0;
%let var=%sysfunc(ifc(&val=1,one,zero));
%put &var;

Una línea, sí, pero una línea muy críptica. ¡Aunque para gustos están los colores!

Nada si tienes los contactos necesarios o te manejas por la parte turbia de internet. Y no tienes inconveniente en mantenerte del lado equivocado de la ley.

O unos 120 dólares si te conformas con una licencia restrictiva de SAS adquiriendo llamada SAS Learning Edition, que puede descargarse desde los servidores de SAS o adquirir en lugares tipo Amazon. Claro, si tus conjuntos de datos no tienen más de 1.500 filas (a fecha de hoy).

Existe un producto específico de SAS, SAS Analytics Pro, que incluye los elemntos más habituales del universo SAS (SAS/BASE, SAS/STAT, SAS/GRAPH, etc.) que sólo cuesta

Tengo entendido que el producto tiene algunas restricciones. Por ejemplo, que sólo puedes utlizarlo para analizar tus propios datos. Si un señor viene, te pide que le calcules la media de unos números, utilizas SAS Analytics Pro para ello y le cobras, estás violando los términos de uso. ¡Cosas veredes!

Pero, ¿cuánto valen los productos de verdad? La información no es pública y está sujeta a variaciones según múltiples criterios: país, tipo de cliente, habilidad del comercial, etc. Varía también según el sistema operativo utilizado, el número de CPU’s, el de usuarios,… Es más opaco, en definitiva, que la telefonía móvil.

Pero hay un cliente muy especial que exige datos públicos: el gobierno de los EE.UU. Éste obliga a todos sus proveedores a revelar sus precios por motivos muy democráticos que no imitamos en España, claro. Así, en esta página puede descargarse este documento que detalla los precios de todos los productos de SAS sólo para el gobierno estadounidense.

No sé si parecen caros o baratos a mis lectores. Pero, por ejemplo, al organismo que gestiona la Seguridad Social en Estados Unidos le parecen tan desaforados que está tratando de encontrar una alternativa a SAS de una manera muy peculiar: busca un proveedor (que no sea el mismo SAS Institute) que le proporcione un intérprete de código SAS.

Seguro que Raúl tiene comentarios la mar de jugosos que hacer a esta entrada.

Hoy aprovecho que pasan dos pájaros por el cielo para pegar un tiro que, seguro, es del interés de mis lectores: voy a utilizar un modelo lineal mixto para estudiar los factores que afectan al rendimiento de una familia de queries de SQL complejas.

El objetivo final es contar con criterios empíricos para la optimización de ciertas queries (siento decir optimización de queries: me obliga a ello la voluntad de que los buscadores me indexen donde más búsquedas se vayan a realizar; por una vez, renegaré del talibán ortográfico que llevo dentro) e, indirectamente, ilustrar con datos distintos de los habituales esta técnica estadística.

La query tiene este aspecto:

select * from ( un carajal de tablas y subconsultas ) where
fecha = :fecha: and
unidad = :unidad: and
cuenta = :cuenta:
;

Los dos factores que se consideran potencialmente críticos para el rendimiento de la query son el número de filas correspondientes a una unidad dada (unit.size en lo sucesivo) en una de las tablas subyacentes y el número de filas correspondientes a una cuenta determinada en otra de dichas tablas (account.size).

Además, se observa que la carga de la base de datos varía mucho en función de factores externos (¿otros usuarios accediendo a ella?) fuera de nuestro control.

Por tanto, analizaremos los tiempos de ejecución de la query en función de unit.size, account.size y, finalmente, de los factores no controlados. Para ello, se selecionamos 200 queries, es decir, 200 combinaciones de unidad, cuenta y fecha. Cada una de las queries se ejecuta 5 veces.

n <- 20
unidades <- sample( unidades, n )
epigrafes <- sample( epigrafes, n )
fechas <- sample( fechas, n, replace = TRUE )

Además, se aleatoriza con respecto al tiempo. Esto se hace para evitar que el efecto de las fluctuaciones de carga del servidor se confundan con el efecto de las distintas queries. Distribuyendo aleatoriamente en el tiempo el momento de la ejecución de una misma query se reduce la posibilidad de que todas las iteraciones de una de ellas se realicen en periodos de carga anormalmente alta o baja (los puristas del diseño experimental nos aplaudirían en este punto sólo a medias).

n.rep <- 5
iter.order <- sample( rep( 1:n, n.rep ) )

Finalmente, se construye el conjunto de datos mediante

salida <- do.call( rbind, sapply( iter.order, foo.sample, simplify = F ) )

donde foo.sample es la función que ejecuta la query contra la base de datos, cuenta el número de registros, etc.

Los datos resultantes pueden descargarse aquí. Se trata de una tabla con cuatro columnas: el identificador de la query, el tiempo de ejecución, unit.size y account.size. (Debería ser ocioso decir aquí que el identificador está asociado al enunciado de la query y no su ejecución puesto que cada una de ellas se repite 5 veces).

Análisis estadístico

Descargamos y normalizamos en primer lugar los datos:

dat <- read.table( url(
   "http://www.datanalytics.com/uploads/query_time_analysis.csv" ),
   header = T ) 

normalize <- function( x ) ( x - mean( x ) ) / sd( x )

dat$unit.size <- normalize( dat$unit.size )
dat$account.size <- normalize( dat$account.size )
dat$id <- factor( dat$id )

(La normalización se hace especialmente para facilitar la interpretación del modelo que se plantea más abajo). A continuación, creamos un objeto de la clase groupedData,

library( nlme )
dat <- groupedData( query.time ~ unit.size + account.size  | id, data = dat )

que viene a ser un data.frame con información sobre cómo ciertas filas están asociadas entre sí. De forma que si uno hace

plot( groupedData( query.time ~ 1 | id, data = dat ) )

se obtiene el siguiente gráfico:

Se aprecia en él cómo la varianza de los tiempos de ejecución crece con éstos. Además, por consideraciones relativas a la construcción de los datos —un cruce de varias tablas de cada una de las cuales se extrae un número variable de filas— hay razones para intuir una estructura multiplicativa en los datos. Eso nos hace considerar el uso de logaritmos. De hecho,

plot( groupedData( log(query.time) ~ 1 | id, data = dat ) )

que produce la gráfica

que tiene mejor aspecto. No es todo lo bueno que uno quisiera, pero tiene mejor aspecto. Nótese además, cómo el que los datos aparezcan ordenados en la figura puede hacer sobreestimar el efecto de la dispersión de la varianza. En realidad, más adelante, se plantea como ejercicio verificar de una manera más canónica cómo el tomar logaritmos no deja de tener sentido.

Finalmente, planteamos el modelo mixto usando la función lme del paquete nlme:

modelo <- lme( log2( query.time ) ~ unit.size + account.size,
   random = ~1 | id , data = dat )

El modelo consta términos fijos (unit.size y account.size) y de una parte aleatoria, ~1 | id. La parte aleatoria, de acuerdo con el consejo recurrente de mi colega Oliver Núñez en r-help-es, comprende aquellos términos del modelo que variarían de realizarse de nuevo el experimento. Y, efectivamente, la próxima vez que se ejecute una query, es improbable que ésta sea una de las 200 seleccionadas más arriba. El enunciado particular de la query es, por lo tanto, variable (o aleatorio, en nuestro contexto).

Alternativamente, usando modelos no mixtos, podría plantearse el modelo equivalente

modelo.lm <- lm( log2( query.time ) ~
   id + unit.size + account.size, data = dat )

que haría aparecer 200 (técnicamente, 199 porque no se ha eliminado el término independiente del modelo) coeficientes nuevos que representarían la variación en los tiempos de ejecución atribuida a cada enunciado de query en particular (variación que puede deberse a la distinta distribución de la filas correspondientes a distintas unidades o cuentas en el disco, etc.). Obviamente, estos coeficientes no son de mayor interés en sí mismos. A lo más, interesa saber si existen variaciones sustanciales entre las distintas queries que pudieran ser indicativas de algún fenómenono tenido en cuenta.

Mediante

summary( modelo )

se obtiene:

Linear mixed-effects model fit by REML
Data: dat
AIC    BIC  logLik
3259.9 3284.4 -1625.0

Random effects:
Formula: ~1 | id
(Intercept) Residual
StdDev:  0.00029632   1.2220

Fixed effects: log2(query.time) ~ unit.size + account.size
Value Std.Error  DF t-value p-value
(Intercept)  0.70249  0.038644 800 18.1785  0.0000
unit.size    0.19907  0.038665 197  5.1486  0.0000
account.size 0.11006  0.038665 197  2.8465  0.0049
Correlation:
(Intr) unt.sz
unit.size    0.000
account.size 0.000  0.008

Standardized Within-Group Residuals:
Min       Q1      Med       Q3      Max
-1.82818 -0.82791 -0.23734  0.70098  3.63443

Number of Observations: 1000
Number of Groups: 200

De la salida anterior interesan varios valores:

• En la sección correspondiente a los efectos aleatorios, la desviación estándar del término independiente, 0.00029632, que parece indicar que apenas hay variación entre las distintas queries.
• En la misma sección, el valor relativamente elevado, 1.2220, del residuo. Eso indica que existe una variación importante entre ejecuciones distintas de la misma query debidas, probablemente, a las distintas condiciones de carga del servidor en el momento de la ejecución. ¡Muy sintomático!
• En la sección de los términos fijos, el término independiente, de valor 0.7 y altamente significativo, debido a la peculiar normalización de los datos, indicaba un tiempo medio de ejecución de 2^0.7 = 1.62 segundos.
• En dicha sección los coeficientes de unit.size y de account.size, también altamente significativos, auguraban un rendimiento desigual de las consultas, posiblemente inaceptablemente desigual.

Finalmente, puede verse un gráfico de diagnóstico del modelo haciendo

plot( modelo )

que produce

y que muestra cómo la varianza de los residuos no parece variar apreciablemente con el tamaño predicho del tiempo de ejecución.

Colofón

En primer lugar, quiero dejar planteados varios ejercicios que serán sin duda de sumo provecho para los más inquietos e interesados de mis lectores:

• Generar y analizar el gráfico de diagnóstico para el modelo análogo en el que no se toman logaritmos de la variable objetivo.
• Calcular la varianza de los coeficientes asociados a cada una de las queries del modelo no mixto modelo.lm y compararla con la del término independiente de la parte aleatoria del modelo mixto.
• Reconstruir el modelo sin normalizar las variables previamente y analizar (comparativamente) los resultados.

Finalmente, doy respuesta a lo que más de uno se estará preguntando: todo esto, ¿para qué? Puede parecer, y alguno así me lo ha manifestado, un ejercicio ocioso. Pero, la verdad, los números muestran varios indicios fundamentales:

• El primero, que parece necesario identificar ventanas temporales en las que el servidor esté desocupado para realizar pruebas de rendimiento: dos desviaciones estándar de ruido intra-query multiplican los tiempos de ejecución en un factor de 5.4264 (=2^(2*1.22) ). Es algo que se intuía pero no se cuantificaba. Y que ponía en entredicho pruebas de rendimiento realizadas anteriormente.
• El segundo, que las dependencias del tiempo de ejecución con respecto al tamaño de las subqueries no son, como se suponía un tanto cándidamente, O(1). Más bien, son O(n*m). Y esto sugirió alterar el orden de los cruces de las tablas para lograr un plan de ejecución alternativo más rápido y robusto (que, de hecho, se encontró).

Desde que dejé de ser uno de ellos, a esa gente que vive en un mundo en el que las cifras tienen un cero de más sólo me la tropiezo en los ascensores. Los oigo hablar de potencias de motores, de la piscina del chalé y de lo mal que está el servicio. Si de verdad tuviesen interés en aquello por lo que les pagan, seguro, leerían esta entrada y no se perderían ni una coma de lo que sigue a continuación.

He de reconocer que no es mío: lo traduzco más o menos libremente de aquí. Está está dirigido a los CIOs de empresas intensivas en SAS y escrito al hilo de la reciente sentencia en el caso que enfrenta a SAS y WPS. Pero quienes toman cortados conmigo por la mañana me han oído alguna vez enumerar las siguientes buenas razones:

1. Dado que SAS se compra/alquila “a la carta”, catalogar todas las licencias de productos de SAS e identificar las aplicaciones que utilizan cada una de ellas. Es increíble cómo muchas compañías —dice el autor y corroboro por experiencia personal— ignoran que han estado pagando por productos de SAS que llevan años sin utilizar.
2. Verificar que los productos de SAS corren sobre las plataformas más económicas, dado que los precios de las licencias de SAS varían sustancialmente en función del hardware y sistema operativo subyacente.
3. Identificar las aplicaciones de SAS más susceptibles de ser migradas a WPS. Dicha migración no es enteramente automática, pero debería ser sencilla en algunos casos, como el de aplicaciones que realizan movimientos de datos mediante pasos data o crean informes. Después de eso, realizar pruebas de concepto sobre dichas aplicaciones en WPS. Y, por supuesto, tener al corriente al comercial de SAS de dichas iniciativas.
4. Probar herramientas modernas de ETL como puedan serlo Pentaho DI (antiguamente Kettle). Plantear tambiéne el uso de otras herramientas de BI open source sólidas como Jaspersoft o Pentaho para la creación de informes, OLAP y cuadros de mando.
5. Identificar nuevos proyectos en los que poder comenzar a probar la capacidad gráfica y analítica de R. R es muy distinto a SAS, por lo que la organización tendría que incorporar consultores —yo me postulo— que pudieran poner en marcha la iniciativa. Pero la satisfacción de los usuarios estaría más que garantizada.

Hace poco, IDC —una empresa que hace estudios de mercado a nivel global de distintas herramientas de sofware y hardware — hizo público su informe periódico Worldwide Business Intelligence Tools 2009 Vendor Shares. En su página 8, la más jugosa del informe, aparece la tabla que reproduzco a continuación:

Puede apreciarse cómo en el segmento de la minería de datos (que viene a ser a lo que se refieren con lo de advanced analytics) es SAS el claro dominador con IBM/SPSS en una débil segunda posición.

Pero la gráfica anterior se refiere, únicamente, a facturación. Y las herramientas de software libre, por muchos usuarios que tengan, nunca van a aparecer en la tabla anterior. Por eso, invito a los lectores de este blog a echarle un vistazo a esta gráfica y a esta otra que he obtenido de este hilo de mensajes de la lista de correo de R.

¿A que se aprecia una tendencia bien distinta?

Ayer quise publicar un comentario a este artículo sobre el reciente veredicto del caso SAS vs. WPS. Aunque ya hablé de eso el otro día, como no me lo publicaron (por ingnotas razones), hago constancia de mi comentario aquí:

The article is misleading in the sense that computer languages (and there exists such a thing as a SAS computer language) are not subject to copyright.

I am free to write my own interpreter/compiler of Java, C, Python, PL/SQL or Logo. And many companies do: Microsoft has implemented their own interpreter of Python, IronPython, etc. Anybody is entitled to write yet another interpreter of SAS code. And, in fact, SAS did not sue WPS for that reason.

Whether the manuals are considered to be so similar that infringe SAS’s copyright can be easily fixed. In fact, you can guide your steps in WPS with widely available SAS’s manuals and books.

Mind that the main allegation of SAS against WPS is a breach of EULA. WPS seemingly acquired a student version of SAS with a restrictive EULA. However, WPS used it for other purposes. However, the court decided that those EULA restrictions were void under European Law: in Europe, law is above private contracts and contracts don’t adhering to the prevailing law are void. For instance, in Europe, I can purchase any piece of software and use it in any legal way that I see fit. I can also buy a car and drive it where I see fit as long as I respect the traffic rules; Fiat, let’s say, would not be able to sell me a car and prohibit me traveling to Portugal, say; or dismantle it and check how the engine works.

WPS acted according to the European law and this is precisely what the courts stated. It is SAS who tries to bind its customers through illegal EULA restrictions.

If the author of the article thinks that it is so easy making a successful clone of a well established software product in Europe and get rich doing it… well, he is quite invited to come here and start cloning!

Resulta curioso cómo ambas partes, WPS y SAS, se arrogan la victoria en el caso. Lo que hace falta ahora es que se abra el mercado y el precio del software se racionalice (¿cero?).