miércoles, 4 de septiembre de 2013

Usando R

Operaciones matriciales usando R

Operaciones matriciales usando R

El paquete Matrix de R contiene 70 clases de matrices y una serie de métodos que nos permiten operar con este tipo de objeto. Recordemos que una matriz es un arreglo, generalmente rectangular, de mxn elementos de un cuerpo K, dispuesto sobre m filas y n columnas. En este caso decimos que la matriz es de orden mxn.

Operaciones con vectores

Un vector, es un arreglo unidimensional de compuesto de n números
v1 <- as.vector(seq(1:10))
v1
##  [1]  1  2  3  4  5  6  7  8  9 10

v2 <- as.vector(seq(11, 20))
v2
##  [1] 11 12 13 14 15 16 17 18 19 20
Podemos operar con vectores
# suma por una constante
v1 + 4
##  [1]  5  6  7  8  9 10 11 12 13 14

# de dos vectores
v1 + v2
##  [1] 12 14 16 18 20 22 24 26 28 30
Los vectores en R no tienen dimensiones, pero sí extensión. No obstante, el transpuesto de un vector genera un vector bidimensional, finalmente, podemos obtener información sobre el contenido del vector, esto el rango del contenido del vector:
dim(v1)
## NULL

length(v1)
## [1] 10

dim(t(v1))
## [1]  1 10

range(v1)
## [1]  1 10
Existen tres posibles maneras de multiplicar vectores en R.
# multiplicacion simple. Los vectores deben tener la misma extensión
v1 * v2
##  [1]  11  24  39  56  75  96 119 144 171 200

v1 * 4
##  [1]  4  8 12 16 20 24 28 32 36 40

# producto externo de un vector para lo que debemos usar el operando %*%t(v)

v1 %*% t(v2)  # devuelve una matriz
##       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
##  [1,]   11   12   13   14   15   16   17   18   19    20
##  [2,]   22   24   26   28   30   32   34   36   38    40
##  [3,]   33   36   39   42   45   48   51   54   57    60
##  [4,]   44   48   52   56   60   64   68   72   76    80
##  [5,]   55   60   65   70   75   80   85   90   95   100
##  [6,]   66   72   78   84   90   96  102  108  114   120
##  [7,]   77   84   91   98  105  112  119  126  133   140
##  [8,]   88   96  104  112  120  128  136  144  152   160
##  [9,]   99  108  117  126  135  144  153  162  171   180
## [10,]  110  120  130  140  150  160  170  180  190   200

# para observar mejor el resultado que se obtiene con esta operacion

v1 %*% t(v1)  #cada elemento de v1 es multiplicado por el elemento correspondiente en v2 con lo que obtenemos las tablas de multiplicar del 1 a l 10
##       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
##  [1,]    1    2    3    4    5    6    7    8    9    10
##  [2,]    2    4    6    8   10   12   14   16   18    20
##  [3,]    3    6    9   12   15   18   21   24   27    30
##  [4,]    4    8   12   16   20   24   28   32   36    40
##  [5,]    5   10   15   20   25   30   35   40   45    50
##  [6,]    6   12   18   24   30   36   42   48   54    60
##  [7,]    7   14   21   28   35   42   49   56   63    70
##  [8,]    8   16   24   32   40   48   56   64   72    80
##  [9,]    9   18   27   36   45   54   63   72   81    90
## [10,]   10   20   30   40   50   60   70   80   90   100
Para extraer, esto es indizar, informacion de un vector, podemos emplear la siguiente notacion v1[1], que nos devuelve el primer elemento; o si queremos extraer el quinto y sexto elemento v1[5:6]
v1[1]
## [1] 1
v1[5:6]
## [1] 5 6
La tercer forma de multiplicar vectores en R multiplica de cada uno de los elementos del v1, en este caso, por los del vector v2 y suma el resultado parcial
v1 %*% v2
##      [,1]
## [1,]  935

Aplicaciones

las operaciones con vectores ayudan a determinar una serie de funciones, que si bien están ya definidas en R, resulta interesante conocer la forma en que se construyen. Por ejemplo la media

# funcion estandar en R, mean()

mean(v1)
## [1] 5.5

# que opera de la siguiente forma: crea vectores de 1 obtiene el producto
# interno de la multiplicacion de los vectores de unos por el vector, en
# este caso, v1 divide el producto interno entre la extensión de los
# vectores de 1

unos <- as.vector(rep(1, length(v1)))  #obtenemos:
unos
##  [1] 1 1 1 1 1 1 1 1 1 1
pInterno <- t(unos) %*% v1  #obtenemos:
pInterno
##      [,1]
## [1,]   55

# media
Media <- pInterno * (1/length(v1))  # y obtenemos el mismo resultado:
Media
##      [,1]
## [1,]  5.5
También se puede llegar a la media por otros caminos:
Media <- t(unos) %*% v1 * (1/length(v1))
Media
##      [,1]
## [1,]  5.5

Media <- t(unos) %*% v1/length(v1)
Media
##      [,1]
## [1,]  5.5

# o simplemente
Media <- sum(v1)/length(v1)
Media
## [1] 5.5
Los vectores también se emplean para calcular la varianza

# la funcion de R
var(v1)
## [1] 9.167
# restamos cada elemento del v1 de la media de v1
v1 - Media
##  [1] -4.5 -3.5 -2.5 -1.5 -0.5  0.5  1.5  2.5  3.5  4.5

# el producto interno de la suma de los cuadrados de la diferencia anterior
# nos da la varianza

varianza <- t(v1 - Media) %*% (v1 - Media) * (1/(length(v1) - 1))
varianza
##       [,1]
## [1,] 9.167
Podemos crear matrices combinando vectores, bien uniendo sus columnas, o bien uniendo sus filas:
# uniendo las columnas cbind()
cbind(v1, v2)
##       v1 v2
##  [1,]  1 11
##  [2,]  2 12
##  [3,]  3 13
##  [4,]  4 14
##  [5,]  5 15
##  [6,]  6 16
##  [7,]  7 17
##  [8,]  8 18
##  [9,]  9 19
## [10,] 10 20

# uniendo las filas rbind()
rbind(v1, v2)
##    [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## v1    1    2    3    4    5    6    7    8    9    10
## v2   11   12   13   14   15   16   17   18   19    20

Vector de datos

En la planificación, operamos con vectores cuando hacemos análisis univariados. Así, por ejemplo, si tenemos datos sobre el cancer de mama en el mundo,dependiendo de la manera en que lo hayamos almacenado, tendremos o un vector columna,nx1 , o un vector fila, 1xn:
a  #vector columna, parte de la informacion sobre cáncer cervical
##    cervical
## 1       2.6
## 2       1.5
## 3       3.4
## 4       0.8
## 5      12.5
## 6       5.4
## 7       0.0
## 8       5.2
## 9       5.5
## 10      1.0
b  #vector fila
##  [1]  2.6  1.5  3.4  0.8 12.5  5.4  0.0  5.2  5.5  1.0
Entre los análisis que normalmente hacemos con vectores, filas, o columnas, están:
Medidas de tendencia central, junto con cuartiles, el valor mínimo y el valor máximo de la distribución:
summary(cervical)  # medidas de tendencia central y cuartiles
##     cervical    
##  Min.   : 0.00  
##  1st Qu.: 2.20  
##  Median : 5.20  
##  Mean   : 6.67  
##  3rd Qu.: 9.50  
##  Max.   :23.70
Observamos que el cancer cervical se distribuye en forma asimétrica, asimetría positiva, con predominio de valores bajos, en 75% de los casos, es decir, de los países en el mundo, presentan tasas de 9.5 por cada cien mil habitantes. Hay, sin embargo, países que presentan tasas inusualmente altas de esta tipología de cancer, hasta 140 veces más altas que el 75% de los casos.
Las medidas anteriores pueden obtenerse en forma individual también:
min(cervical)  # valor mínimo
## [1] 0
quantile(cervical[, 1], probs = 0.25)
## 25% 
## 2.2
median(cervical[, 1])
## [1] 5.2
mean(cervical[, 1])  # media
## [1] 6.671
quantile(cervical[, 1], probs = 0.75)
## 75% 
## 9.5
max(cervical)
## [1] 23.7

Análisis de dispersión
rango <- max(cervical) - min(cervical)  # rango(o recorrido) de la distribución
rango
## [1] 23.7
IQR(cervical[, 1])  # desviación inter-cuartílica
## [1] 7.3
# la varianza y la desviación estandard tienen poco sentido para esta
# distribución del cancer cervical porque esta districión es asimétrica
# (vale decir lo mismo para la media). Su cálculo se puede obtener, en todo
# caso usando las funciones var() y sd()

var(cervical[, 1])
## [1] 30.94
sd(cervical[, 1])
## [1] 5.562

Visualizaciones

El gráfico de caja y bigotes nos ofrece información visual sobre los valores mínimos, máximos y cuartiles de la distribución:
boxplot(cervical, ylab = "Casos por cada cien mil habitantes", main = "Tasas de cáncer cervical en el mundo")
grid()
plot of chunk bp
Puesto que el 75% de los casos presentan valores bajos de incidencia de cáncer de cervical (aquí bajo quiere decir en comparación con todos los valores de esta distribución, y no intenta ofrecer alguna evalución con respecto al cáncer cervical como tal), por ello la media se sesga hacia estos valores bajos. El gráfico de caja y bigote nos indica además que, efectivamente, los valores del cáncer cervical se concentran hacia los valores bajos, de hecho la mediana se ubica en 5 casos por cada cien mil habitantes, en tanto que hay una mayor dispersión hacia los valores altos. Cinco paises presentan tasas de incidencia de cáncer cervical superior a 20 casos por cada cien mil habitantes.
Un gráfico bastante sencillo pero que nos ayuda a comprender rápidamente la distribución de nuestro vector es el grafico de tallo y hoja:
stem(cervical[, 1])
## 
##   The decimal point is at the |
## 
##    0 | 0045555778899
##    1 | 00000011112222333344445555666889
##    2 | 01133344455667999
##    3 | 011124444558
##    4 | 0000133455667799
##    5 | 00112245556799
##    6 | 00123455667
##    7 | 0234667779
##    8 | 0011113334789
##    9 | 0033555678999
##   10 | 1459
##   11 | 1235788
##   12 | 5589
##   13 | 3
##   14 | 01366789
##   15 | 12
##   16 | 3
##   17 | 8
##   18 | 9
##   19 | 77
##   20 | 3
##   21 | 014478
##   22 | 4
##   23 | 7
El gráfico nos confirma la información del gráfico de caja y bigotes. La distribución del cáncer cervical se concentra en los valores bajos, específicamente entre los 0 y 9 casos por cada cien mil habitantes. Adicionalmente la data parece ser por lo menos bimodal. Otro gráfico que podemos elaborar para comprender la data, es un gráfico de barras
barplot(cervical[, 1], ylim = c(0, 25), border = NA)
grid()
plot of chunk unnamed-chunk-2
El gráfico de barras nos sugiere que
  • existen conglomerados de países según la incidencia del cáncer cervical, tal como veíamos con el gráfico de tallo y hoja.
  • posiblemente podemos clasificar a estos países en por lo menos tres grupos: incidencia alta, incidencia media e incidencia baja.
Ordenando la data de menor a mayor podemos apreciar mejor estos posibles conglomerados:
barplot(sort(cervical[, 1]), ylim = c(0, 25), border = NA)
grid()
plot of chunk unnamed-chunk-3
Para verificar la multimodalidad de la distribución, creamos un gráfico de densidad:
plot(density(cervical[, 1], , kernel = "epanechnikov"))
grid()
plot of chunk densidad
La data está sesgada a la derecha, con una concentración hacia los valores bajos, como hemos visto, y se pueden detectar tres conglomerados de a partir de las tasas de cáncer cervical. tenemos por lo menos tres grupos de países de acuerdo con la tasa de cáncer cervical. El análisis del vector de datos sobre cáncer cervical a nivel mundial nos indica que, en general:
  • incidencia de cáncer cervical, en el 75% de los casos, es de 9.5 por cada cien mil personas o menos,
  • podemos detectar por lo menos tres grupos de países de acuerdo con la tasa de incidencia del cancer cervical
  • hay un grupo de países, pequeño, que presentan una incidencia del cáncer cervical 140% más grandes con respecto al 75% de los casos, en los que la incidencia es de 9.5 casos por cada cien mil habitantes o menos

Dos vectores de datos

Queremos analizar en conjunto no sólo la incidencia de cáncer cervical, sino también los países en donde ocurre esto. Para ello vamos incluir el vector columna correspondiente a los países en donde ocurre el hecho:
head(Cervical)
##               pais cervica
## 1      afghanistan     2.6
## 2          albania     1.5
## 3          algeria     3.4
## 4          andorra     0.8
## 5           angola    12.5
## 6 antigua and bar.     5.4
dim(Cervical)
## [1] 191   2
Tenemos entonces 2 columnas, país y cervica, y 191 observaciones. Visualicemos esta relación
plot(sort(Cervical[, 2]), ylab = "tasa por cada cien mil habitantes")
grid()
plot of chunk unnamed-chunk-4
100 paises presentan una incidencia de cáncer de cervical de 5 casos por cada cien mil habitantes. Y 150 países presentan tasas de hasta 10 casos por cada cien mil habitantes, en tanto que aproximadamente 20 países presentan tasas de incidencia superiores a 15 casos por cada cien mil habitantes. Añadamos una tercera columna en la que clasificamos a los países en algunas regiones
regiones[1:10]
##  [1] menaga    europe    menaga    europe    ssa       caribbean menaga   
##  [8] sa        europe    aao      
## Levels: aao asia ca caribbean europe menaga na sa ssa
Observemos cuantos países por regiones tenemos en la data:
reg <- table(regiones)
reg
## regiones
##       aao      asia        ca caribbean    europe    menaga        na 
##        16        27         7        14        44        23         2 
##        sa       ssa 
##        12        46
# agregamos ahora la columna region a nuestrdos dos vectores

regionCervical <- cbind(Cervical, regiones)
head(regionCervical)
##               pais cervica  regiones
## 1      afghanistan     2.6    menaga
## 2          albania     1.5    europe
## 3          algeria     3.4    menaga
## 4          andorra     0.8    europe
## 5           angola    12.5       ssa
## 6 antigua and bar.     5.4 caribbean
barplot(reg, border = NA, las = 2, main = "Regiones en la data", ylim = c(0, 
    50))
grid()
plot of chunk unnamed-chunk-5
Visualizamos ahora la incidencia de cáncer cervica por región
boxplot(regionCervical[, 2] ~ regionCervical[, 3], las = 2, ylim = c(0, 25))
grid()
plot of chunk unnamed-chunk-6

No hay comentarios: