#!/usr/bin/perl

use 5.043.6;
use utf8;

my $version = 0.06;    # as of jue 25 dic 2025 11:18:22 CET

my $t0 = time();

# NOTA de versión: hay dos cosas que no tuve en cuenta hasta
# ahora 1) que la rareza se puede emplear como una dimensión más
# o etiqueta del dato 2) el elefante en la cacharrería es que 
# la etiqueta del dato puede ser el tiempo y mandar en todo


################################################################




# Vamos a suponer que tenemos N datos de M variables
# no ordenados por tiempo ni por nada, sólo por el conjunto
# de datos de una observación

# PASO 0: descriptivas de datos

# PASO 1

# Un primer paso es la caracterización de la observación
# con respecto a las medias individuales por variable, 
# la distancia multidimensional. Se puede prescindir de 
# la raíz

# PASO 2

# Dependiendo de N, también puede ser interesante saber
# en que 'decil' esta cada variable. Esto es un número.
# Si N es muy grande, pueden ser partes por millón. (*)

# Como el número de dimensiones es supuestamente alto,
# hay que prescindir por el momento de supuestos de 
# normalidad y meterlo sólo si es necesario por algo o 
# para algo

# PASO 3

# Que el máximo y mínimo de cada variable es interesante 
# parece lógico, por la que ya tenemos 2*M observaciones 
# mínimas interesantes. Puede haber empates y aumentar 
# el número 

# Esto lleva a que puede ser que interese el número de 
# rango en un orden u otro. Hay que elegir algo consistente,
# por lo que puede ser de más raro a menos o al revés.
# Decidir esto no es sencillo. El objetivo es sacar datos
# interesantes de datos no interesantes

# PASO 4

# El supuesto más general es que los datos son continuos,
# por lo que el valor más frecuente puede no tener mucho
# sentido, salvo en el intervalo medido (*) 

# vie 19 dic 2025 03:42:32 CET

# Cosa suya es reiterar que los datos no tienen marca de
# tiempo ni se trata de una ventana de datos, sino de una
# muestra desreferenciada, y se estudia algo desconocido.

# La semimatriz de covarianzas M*M puede aportar algo, pero
# antes hay que exprimir la relación entre M y N. Cuántos
# datos N se pueden tener (no todos válidos, p.ej.) por cada
# variable medida de M que puede ser 0% o 100% decisiva.
# Parece que va a haber que hacer supuestos sobre lo ya
# decidido (*)

# PASO 5

# Parece interesante relacionar máx-min con N. El oden 
# de magnitud de la diferencia y el orden de magnitud 
# de los datos. Puede haber variables abiertas o cerradas,
# por así decirlo. El valor crítico parece ser log neperiano
# de (max-min) y para ponerlo en relación, /N, o (N-1), en
# purismo

# El caso es que va a haber variables mejor medidas que otras,
# con la muestra N. Como no se ha supuesto ningún ranking de
# preferencias o interés, el interés debe ser inverso a la 
# buena medición  

# No veo sentido a reescalar los datos de ninguna manera, a no
# ser en congujación con lo indicado por aparatos de medida u 
# otro tipo de error sistemático observado y a corregir de la
# manera que sea (redondeo, truncamiento, etc.)

# PASO 6 

# Las distancias son 1D, luego se pueden analizar como una variable
# más. A ver cómo las meto.  

# pequeños datos iniciales para probar
my $N = 100_000;
say "Número de datos = $N";
say "Número de variables = 1000";
say "Número de defectos de datos = 200";
my @data;
for my $linea (0 .. $N){
    $data[$linea][$_] = $_ * rand() * 100 for (1 .. 1000);
    $data[$linea][int(100*rand())] = undef for (1 .. 200);  # datos malos
}

# say "Paso 0\n";
# say "$d datos leídos";

# say "Paso 1\n";
say "\nLista de distancias multidimensionales por línea";

my @mu;
for my $j ( 0 .. scalar(@{$data[0]})-1){
    $mu[$j] = 0;
}

my $d = 0;
my @M;
my @names;
my @nodata;
while ( $d <= $N ){
    for my $j (0 .. scalar(@{$data[0]})-1){
	$nodata[$j] = 0;
	if ( defined $data[$d][$j] ) {
	    $mu[$j] += $data[$d][$j] ;    
	}else{
	    $nodata[$j]++;         # esto se usa después (o no :-)
	}
	$names[$j] = "X" . $j;
	$M[$j] = scalar( @{$data[0]} ) - $nodata[$j] -1;
    }
    $d++;
}
$d--;

say "\na) Lista de medias";
for my $j (0 .. scalar(@{ $data[0] })-1){
    $mu[$j] /= $M[$j];
    say "$names[$j] media $mu[$j]";
}

say "\nb) Lista de distancias cuadráticas a las medias";
my @dist;
for my $i (0 .. $d){
    for my $j (0 .. scalar(@{$data[$i]})-1){
	if (defined $data[$i][$j]){
	    $dist[$i] += ($data[$i][$j]-$mu[$j])**2;
	}
    }
    $dist[$i] = sqrt($dist[$i]);  # esto se mete pq el número puede ser GRANDE
    say "dato $i) dist_m $dist[$i]";          # no porque afecte a la ordenación de distancias
}

# Paso 2
# Lo de los centiles de variables

# say "Paso 2\n";
say "\nCentiles de variables: número de centil por dato";
say "y frecuencia del centil de toda la masa muestral";

my (@max, @min);
@max = (-9e99) x scalar (@{ $data[0] });                        # algún dato habrá
@min = (9e99) x scalar (@{ $data[0] });
    
for my $i (0 .. $d){
    for my $j (0 .. scalar(@{ $data[0]})-1){
	if (defined $data[$i][$j]){
	    $max[$j] = ($data[$i][$j] > $max[$j]) ? $data[$i][$j] : $max[$j];
	    $min[$j] = ($data[$i][$j] < $min[$j]) ? $data[$i][$j] : $min[$j];	     
	}
    }
}
	    
# El criterio de ordenación ascenedente o descendente no se toca,
# de momento, porque sino esto no lo entiendo ni yo. Ascendente.

my @centiles;  # va a ser un array [0..100][$j número de variable] 

for my $j (0 .. scalar(@{ $data[0] }) -1 ){
    for my $z (0 .. 100){ 
	$centiles[$z][$j] = $min[$j] + ($max[$j] - $min[$j]) * ($z /100) ;
	# no sé si esto es correcto...
	# el centil 0 sólo contiene el mínimo, parece bien
	# el centil 100 sólo contiene el máximo, parece bien
	# puede haber @frecuencias mayores por repetición
    }
    $centiles[101][$j] = $max[$j];     # un truco para ir tirando, podría ser = $centiles[100][$j] 
}

my @frecuencia;      
                 # aquí hay que decidir las dimensiones, [$i][$i+1] ???
                 # no tomar atajos: [$i centil][$j variable] semimatriz

for my $j ( 0 .. scalar(@{ $data[0] })-1 ){
    for my $z (0 .. 100){
	$frecuencia[$z][$j] = 0;
    } 
    for my $k (0 .. $d){
	if ( defined $data[$k][$j] ){
	    for my $z (0 .. 100) {
		if ( $data[$k][$j] >= $centiles[$z][$j] && $data[$k][$j] < $centiles[$z+1][$j] ){
		    $frecuencia[$z][$j]++;                                         # no acumulada, ya se verá lo acumulado
		}
	    }
	}
    }
}


for my $j (0 .. scalar(@{ $data[0] })-1){
    for my $z (0 .. 100) {
	print "$names[$j] - centil $z% : $centiles[$z][$j] ";
	print "-" x ( $frecuencia[$z][$j]/100 ), "\n";
    }
}

# NOTA IMPORTANTE: elegí c1 > dato <= c2 porque todo es menor que 100% pero poco mayor que 0%
# La cosa es que los semisegmentos deben poder sumarse, y así se hace. 

say "\nLista de número de datos por amplitud rango (log10) inverso";
my @gradocontrol;                             # sobre variables!!!
for my $j ( 0 .. scalar(@{ $data[0] })-1 ){
    if ($max[$j] == $min[$j]){
	say "ERROR: DATO de variable $j ES PLANO" ;
	$gradocontrol[$j] = "log10 = +inf xor -inf, es plano";
	next;
    }
    $gradocontrol[$j] = ( $N-$M[$j] ) / log10 ( abs($max[$j] - $min[$j]) );   # resto M por grados de libertad, pero sale grande 
    # sobre si $max - $min puede dar negativo, como que no
    # sobre si $max - $min puede dar 0 -> excepción dato plano
}
for my $j (0 .. scalar( @{ $data[0] })-1 ) {
    say "$names[$j] -> amplitud $gradocontrol[$j]";
}

# Ahora tengo que hacer subrutinas de lo de arriba para poder darles
# la lista que quiera, que interesa que sea la de las distancias como 
# cualquier otra variable, y de la diferencia de las distancias si no
# fuese porque la muestra de observaciones bien puede estar desordenada

# O no. Se repite y fuera. Parte B es que se podría hacer A. cluster,
# pero eso puede no tener sentido para todo tipo de datos. Más importante
# es poner etiquetas a las variables X0 .. X_N

##################################### SIGIENTE FASE lun 22 dic 2025 23:15:45 CET

my @cluster;                  # NO LO NORMAL, PORQUE NO ES NI 2D
my @frees = 0 .. $N;
my @freesord = sort { $dist[$a] <=> $dist[$b] } @frees;

my @diffree;
$diffree[0] = 0;   # necesario para que el tamaño de las listas sea compatible
for my $i (1 .. scalar(@freesord)-1){
    $diffree[$i] = $dist[$freesord[$i]] - $dist[$freesord[$i-1]];
}

my @diffreeord = sort { $diffree[$a] <=> $diffree[$b] } @frees;

say "Orden de diferencias:";  
for my $i ( @diffreeord ){
    say "orden $i) diferencia $diffree[$i] -> punto original $freesord[$i]";
}


goto planB;

my @tmp;
my @nclus = 0;                    # es una lista de 'flags', podría ser booleana. Sí, mejor
my $nlist = 0;                    # número incremental, siempre hacia arriba, de la lista-resultado 
my $n = 0;                        # el contador de grupos activos
my $iter = 0;

do {
   # ...
   # - equiprobabilidad
   # check min: inside @frees, inside @clusters, @frees to @clusters (1D, no centroids)
    
   # MIN 3 groups: a) easy, min dist
   # b) easy (1D to 1D lists) no hay intermedios
   # c) punto a segmento, mínimo  

   # por el caso b) se necesita una ordenación previa 
    
   # si hay ordenación, tonces directamente se coge la diferencia y se agrupa siempre por mínima diferencia  
    
   # ... 
    
    
   # - con peso de rareza (lo del @gradocontrol)
   # Esto ya lía las cosas

    $iter++;
    
    my ($candidate1, $candidate2, $candidate2b, $candidate3, $candidate3b) = (0,0,0,0,0);  
    my $mini = 9e99;
    for my $f (@frees){
	my $g = $dist[$f];
	if ($g < $mini){
	     $mini = $f;
	}
    }
    $candidate1 = $mini;                 # un punto que pasa a ser origen de un grupo
    
    my $mini1 = 9e99;
    my $mini2  = 9e99;
    if ($n == 0){ goto ya; }
    for my $c1 (1 .. scalar(@cluster)-1) {
	next unless ($nclus[$c1] > 0);
	for my $c2 ($c1 .. scalar(@cluster)-1) {
	    next unless ($nclus[$c2] > 0);
  
	    for my $c ( @{ $cluster[$c2] } ){
		for my $content ( @{ $cluster[$c1] } ){
		    if ($dist[$content] < $mini1 && $dist[$c] < $mini2){
			$mini1 = $c1;        
			$mini2 = $c2;
		    }
		}
	    }
	    
	}	
    }
ya:    
    $candidate2 = $mini1;              # dos grupos que se unen, el número de grupos baja uno             
    $candidate2b = $mini2;

    $mini1 = 9e99;
    $mini2 = 9e99;
    for my $c1 (1 .. scalar(@cluster)-1){
	for my $c ( @{ $cluster[$c1] } ){
	    for my $f (@frees){
		if ($dist[$c] < $mini1 && $f < $mini2){
		    $mini1 = $c1;
		    $mini2 = $f;
		}
	    }
	}
    }
    $candidate3 = $mini1;              # un punto se une a un grupo, no cambia el número de grupos
    $candidate3b = $mini2;
    # new cluster whatever comes
#    $nlist++;                   # 0j0, el primero es 1, no 0
    if ($candidate1 < $candidate3 && $candidate1 < $candidate3b){
	push @{ $cluster[$nlist] }, $candidate1;       # caso de creación
	$nclus[++$nlist] = 1;
	$n++;                                            # se añade 1
    }elsif ($candidate3 < $candidate2 && $candidate3b < $candidate2b){
	#  @nclus igual 
	push @{ $cluster[$candidate3] }, $candidate3b;
    }else{
	@tmp = @{ $candidate2b };
	push @{ $cluster[$candidate2] }, @tmp;                # se quitan los uno viejo               
#	$nclus[$candidate2] = 0;                   #
	$nclus[$candidate2b] = 0;                        # es decir, cada cosa va aparte
#	$nclus[$nlist] = 1;
	$n += - 1;                                     # resumen
    }

    # say blah    -> sólo se lista lo activo y de lo libre, sólo el número 
#    system "clear";
    say "\nIter $iter";
    say "Número de libres = ", scalar(@frees);
    say "Número de grupos activos = $n";
    for my $i (1 .. $nlist){
	if ($nclus[$i] > 0){
	    say "grupo $i) ", @{ $cluster[$i] };
	}
    }
}until (scalar(@frees) == 0 && $n == 1); 

planB:

say "\nStop program $0 at ", scalar(localtime());

say "Elapsed Time ", int( 10 * (0.5 + (time() - $t0) / 60)) /10 , " minutes";

exit 3;


sub log10 {
    my $argu = shift;
    return (log ($argu)) / (log (10));
}
    


__END__
  
Faltan cosillas
v.1.07 - crear una distancia de *rareza* que sea laqa media ponderada
  de los datos por el peso del @gradodecontrol relativo de variable,
  por ejemplo. No obstante, ya se mide con última (max) @dist