e2001 <- dget("/home/etptupaf/datos/e2001.frame")
Vizcaya <- e2001[e2001$Territorio==48,9:18]
Vizcaya[is.na(Vizcaya)] <- 0

dibujaXY <- function(x,y,tipo="no",etiquetas,valores,col=3,pch=".") {
    #
    #  Esta función auxiliar imprime etiquetas, caracteres, círculos
    #  o estrellas en las coordenadas x,y. Si queremos etiquetas, debemos
    #  proporcionar el argumento etiquetas; si queremos circulos o estrellas,
    #  el argumento valores (que es distinto para círculos y estrellas;
    #  véase la documentación de la función symbols).
    #
    switch(tipo,
           texto = text(x,y,labels=etiquetas,col=col),
           corto = text(x,y,labels=substring(etiquetas,1,3),col=col),
           punto = points(x,y,col=col,pch=pch),
           circulo = {valores <- sqrt(valores/pi);
                      symbols(x,y,circles=valores,inches=0.7,add=TRUE)},
           estrella = {n <- length(valores);
                        valores <- matrix(c(valores,0.2*valores),n,10);
                        symbols(x,y,stars=valores,inches=0.7,add=TRUE)},
           no = 
       )
}

"ac" <-
  function(x,filas=c("texto","punto","circulo","estrella","no"),
             columnas=c("texto","punto","circulo","estrella","no"), nv = 2) 
{
  #
  #   Analisis de correspondencias simple.
  #
  #  x        = tabla de contingencia.
  #  filas    = Modo en que queremos representar los puntos fila.            
  #  columnas = Modo en que queremos representar los puntos columna.
  #  nv       = número de ejes que queremos representar. Casi invariablemente
  #             dejaremos también este argumento en paz, y obtendremos una
  #             representación en dos dimensiones con ayuda de los dos
  #             primeros vectores propios. Si fijáramos nv=3, tendríamos las
  #             representaciones (Eje1,Eje2), (Eje1,Eje3), (Eje2,Eje3).
  #
  x <- as.matrix(x) / sum(x)
  #
  # Ahora "x" ha sido reducida a frecuencias relativas: es la F  del
  # desarrollo en los apuntes o realizado en clase.
  #
  marfilas <- apply(x, 1, sum)
  marcolum <- apply(x, 2, sum)
  #
  # Sumamos respecto de filas y columnas. "marfilas" es el vector de
  # frecuencias marginales f(i,.) y "marcolum" el de frecuencias marginales
  # f(.,j).
  #
  descomp <- svd(x/outer(sqrt(marfilas), sqrt(marcolum)))
  #
  # La matriz de la que se calcula la descomposición en valores singulares
  # en la linea precedente es la Z del análisis en clase. La SVD da de una
  # vez (y de forma computacionalmente muy económica) los vectores propios
  # de ZZ' y de Z'Z que necesitamos.
  #
  # Los elementos descomp$u y descomp$v de la descomposición en valores
  # singulares son los vectores propios A y B del desarrollo en clase.
  #
  R <- (x/outer(marfilas, sqrt(marcolum))) %*% 
    (descomp$v %*% diag(descomp$d))[, -1]
  C <- (t(x)/outer(marcolum, sqrt(marfilas))) %*% 
   (descomp$u %*% diag(descomp$d))[, -1]
  #
  #  ¡Ya está! Las lineas que siguen, meramente hacen los gráficos.
  #  Si nv = 1, se hace un sólo gráfico con los ejes 1 y 2.
  #
  for(i in 1:(nv - 1)) {
    for(j in (i + 1):nv) {
      #
      #  Cuando hacemos el gráfico de los Ejes i,j, las coordenadas de
      #  los puntos que representamos salen de las columnas correspondientes
      #  de R y C.
      #
      puntosx <- c(R[, i], C[, i])
      puntosy <- c(R[, j], C[, j])
      #
      #  La línea siguiente emplea la técnica habitual del "gráfico mudo",
      #  en el que fijamos las escalas y rotulamos los ejes para luego añadir
      #  los puntos
      #
      plot(puntosx,puntosy,type="n",
           xlab = paste("Eje ", i), ylab = paste("Eje ", j))
      #
      #  Dependiendo de los valores de "filas" y "columnas"
      #  escribimos en las coordenadas adecuadas las etiquetas o símbolos
      #  de los puntos fila y columna. (Podría interesarnos hacer 
      #  filas="no", si hubiera tantos puntos fila que oscurecieran el
      #  gráfico).
      #
      dibujaXY(C[,i],C[,j],tipo=columnas,etiquetas=dimnames(x)[[2]],
               valores=marcolum,col=2,pch=".")
      dibujaXY(R[,i],R[,j],tipo=filas,etiquetas=dimnames(x)[[1]],
               valores=marfilas,col=3,pch=".")
    }
  }
  #
  #  Despues de hacer los gráficos requeridos, la función todavía retorna
  #  las coordenadas de filas y columnas por si queremos hacer algún otro
  #  análisis o gráfico con ellas.
  #
  return(list(R, C))
}
#
#   Ejemplos de uso
#
a <- ac(Vizcaya,filas="texto",columnas="circulo")
a <- ac(Vizcaya,filas="texto",columnas="estrella")
a <- ac(Vizcaya,filas="corto",columnas="estrella")
a <- ac(Vizcaya,filas="circulo",columnas="estrella")