Regresión lineal en NumPy con matrices muy grandes: ¿cómo ahorrar memoria?

Así que tengo estas matrices gigantescas X e Y. X e Y tienen 100 millones de filas, y X tiene 10 columnas. Estoy tratando de implementar la regresión lineal con estas matrices, y necesito la cantidad (X^T*X)^-1 * X^T * Y ¿Cómo puedo calcular esto de la manera más eficiente en el espacio posible?

Ahora mismo tengo

 X = readMatrix("fileX.txt") Y = readMatrix("fileY.txt") return (X.getT() * X).getI() * X.getT() * Y 

¿Cuántas matrices se almacenan en la memoria aquí? ¿Se almacenan más de dos matrices a la vez? Hay una mejor manera de hacerlo?

Tengo alrededor de 1.5 GB de memoria para este proyecto. Probablemente pueda estirarlo a 2 o 2.5 si cierro todos los otros progtwigs. Idealmente, el proceso también se ejecutaría en un corto período de tiempo, pero el límite de memoria es más estricto.

El otro enfoque que he intentado es guardar los pasos intermedios del cálculo como archivos de texto y volver a cargarlos después de cada paso. Pero eso es muy lento.

el tamaño de X es 100e6 x 10 el tamaño de Y es 100e6 x 1

así que el tamaño final de (X^T*X)^-1 * X^T * Y es 10 x 1

Puedes calcularlo siguiendo el siguiente paso:

  1. calcular a = X^T*X -> 10 x 10
  2. calcule b = X^T*Y -> 10 x 1
  3. calcula a^-1 * b

Las matrices en el paso 3 son muy pequeñas, por lo que solo necesita realizar algunos pasos intermedios para calcular 1 y 2.

Por ejemplo, puede leer la columna 0 de X e Y, y calcularla por numpy.dot(X0, Y) .

para el tipo de dato float64, el tamaño de X0 e Y es de aproximadamente 1600M, si no se ajusta a la memoria, puede llamar numpy.dot dos veces para la primera mitad y la segunda mitad de X0 e Y por separado.

Entonces, para calcular X^T*Y necesita llamar numpy.dot 20 veces, para calcular X^T*X necesita llamar numpy.dot 200 veces.

Una propiedad interesante de la regresión de mínimos cuadrados ordinarios es que si tiene dos conjuntos de datos X1, Y1 y X2, Y2 y ya ha calculado todos

  • X1 ‘* X1
  • X1 ‘* Y1
  • X2 ‘* X2
  • X2 ‘* Y2

Y ahora desea hacer la regresión en el conjunto de datos combinado X = [X1; X2] e Y = [Y1; Y2], en realidad no tienes que volver a calcular mucho. Las relaciones

  • X ‘* X = X1’ * X1 + X2 ‘* X2
  • X ‘* Y = X1’ * Y1 + X2 ‘* Y2

espera, así que con estos computados solo calculas

  • beta = inv (X ‘* X) * (X’ * Y)

y tu estas listo. Esto conduce a un algoritmo simple para OLS en conjuntos de datos muy grandes:

  • Cargue en parte del conjunto de datos (por ejemplo, el primer millón de filas) y calcule X ‘* X y X’ * Y (que son matrices bastante pequeñas) y guárdelas.
  • Siga haciendo esto durante el próximo millón de filas, hasta que haya procesado todo el conjunto de datos.
  • Sume todos los X ‘* Xs y X’ * Ys que haya almacenado
  • Calcular beta = inv (X ‘* X) \ (X’ * Y)

Eso no es apreciablemente más lento que cargar todo el conjunto de datos a la vez y usa mucha menos memoria.

Nota final: nunca debe calcular la versión beta calculando primero (X ‘* X) y encontrando su inverso (por dos razones: 1. es lenta y 2. es propensa a errores numéricos).

En su lugar, deberías resolver el sistema lineal –

  • (X ‘* X) * beta = X’ * Y

En MATLAB se trata de un sencillo de una sola línea.

 beta = (X' * X) \ (X' * Y); 

y espero que numpy tenga una forma similar de resolver sistemas lineales sin necesidad de invertir una matriz.

La memoria RAM es bastante barata, deberías considerar invertir. Un sistema con 24 Gig de RAM no necesariamente cuesta un arm y una pierna más, uno de los servidores de gama baja de Dell puede incluir tanto.

Si las matrices son dispersas (muchos ceros), use una clase de matriz dispersa para ahorrar una gran cantidad de RAM.

Si las matrices no son dispersas, querrá más memoria RAM (o al menos más memoria virtual), o realizar sus operaciones de matriz utilizando archivos de disco.

Los archivos de disco son, por supuesto, un orden de magnitud más lento que la memoria RAM, y el sistema de memoria virtual podría ser peor que eso, según los patrones de acceso.