1 #Author: Pierre Ratinaud
2 #Copyright (c) 2008-2020 Pierre Ratinaud
5 pp<-function(txt,val) {
9 MyChiSq<-function(x,sc,n){
11 E <- outer(sr, sc, "*")/n
12 STAT<-sum(((x - E)^2)/E)
16 MySpeedChi <- function(x,sc) {
18 E <- outer(sr, sc, "*")
19 STAT<-sum((x - E)^2/E)
23 find.max <- function(dtable, chitable, compte, rmax, maxinter, sc, TT) {
24 ln <- which(dtable==1, arr.ind=TRUE)
26 lo[1:nrow(dtable)] <- 0
27 for (k in 1:nrow(ln)) {lo[[ln[k,1]]]<-append(lo[[ln[k,1]]],ln[k,2])}
28 for (k in 1:nrow(dtable)) {lo[[k]] <- lo[[k]][-1]}
29 ## lo<-lo[-c(1,length(lo))]
31 ## compte <- compte + 1
32 ## chitable[1,l]<-chitable[1,l]+1
33 ## chitable[2,l]<-chitable[2,l]-1
34 ## chi<-MyChiSq(chitable,sc,TT)
35 ## if (chi>maxinter) {
42 chi<-MyChiSq(chitable,sc,TT)
48 chitable[1,l]<-chitable[1,l]+1
49 chitable[2,l]<-chitable[2,l]-1
51 res <- list(maxinter=maxinter, rmax=rmax)
59 CHD<-function(data.in, x=9, mode.patate = FALSE, svd.method, libsvdc.path=NULL){
60 # sink('/home/pierre/workspace/iramuteq/dev/findchi2.txt')
62 row.names(dataori) <- rownames(data.in)
64 colnames(dtable) <- 1:ncol(dtable)
67 pp('ncol entree : ',ncol(dtable))
68 pp('nrow entree',nrow(dtable))
72 print('vire colonnes vides en entree')#FIXME : il ne doit pas y avoir de colonnes vides en entree !!
75 dtable<-dtable[,-which(sdt==0)]
76 print('vire lignes vides en entree')
79 rowelim<-as.integer(rownames(dtable)[which(sdt==0)])
80 print('&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&')
82 print('&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&')
83 dtable<-dtable[-which(sdt==0),]
88 listmere[[clnb]]<-mere
89 listmere[[clnb+1]]<-mere
90 list_fille[[mere]] <- c(clnb,clnb+1)
91 listcol[[clnb]]<-vector()
92 listcol[[clnb+1]]<-vector()
93 #extraction du premier facteur de l'afc
95 pp('taille dtable dans boucle (col/row)',c(ncol(dtable),nrow(dtable)))
96 afc<-boostana(dtable, nd=1, svd.method = svd.method, libsvdc.path=libsvdc.path)
97 pp('SV',afc$singular.values)
98 pp('V.P.', afc$eigen.values)
99 coordrow <- as.matrix(afc$row.scores[,1])
100 coordrowori<-coordrow
101 row.names(coordrow)<-rownames(dtable)
102 coordrow <- cbind(coordrow,1:nrow(dtable))
103 print('deb recherche meilleur partition')
104 ordert <- as.matrix(coordrow[order(coordrow[,1]),])
105 ordert <- cbind(ordert, 1:nrow(ordert))
106 ordert <- ordert[order(ordert[,2]),]
110 dtable <- dtable[order(ordert[,3]),]
111 sc <- colSums(dtable)
114 sc2 <- colSums(dtable) - sc1
115 chitable <- rbind(sc1, sc2)
120 inert <- find.max(dtable, chitable, compte, rmax, maxinter, sc, TT)
121 print('@@@@@@@@@@@@@@@@@@@@@@@@@@@@')
122 pp('max inter phase 1', inert$maxinter/TT)#max(listinter))
123 print('@@@@@@@@@@@@@@@@@@@@@@@@@@@@')
124 ordert <- ordert[order(ordert[,3]),]
125 listclasse<-ifelse(coordrowori<=ordert[(inert$rmax),1],clnb,clnb+1)
126 dtable <- dtable[order(ordert[,2]),]
129 #dtable<-cbind(dtable,'cl'= as.vector(cl))
131 N1<-length(listclasse[listclasse==clnb])
132 N2<-length(listclasse[listclasse==clnb+1])
135 ###################################################################
136 # reclassement des individus #
137 ###################################################################
143 ln <- which(dtable==1, arr.ind=TRUE)
145 lnz[1:nrow(dtable)] <- 0
147 for (k in 1:nrow(ln)) {lnz[[ln[k,1]]]<-append(lnz[[ln[k,1]]],ln[k,2])}
148 for (k in 1:nrow(dtable)) {lnz[[k]] <- lnz[[k]][-1]}
151 while (malcl!=0 & N1>=5 & N2>=5) {
153 listsub[[it]]<-vector()
154 txt <- paste('nombre iteration', it)
155 #pp('nombre iteration',it)
158 t1<-dtable[which(cl[,1]==clnb),]#[,-ncol(dtable)]
159 t2<-dtable[which(cl[,1]==clnb+1),]#[,-ncol(dtable)]
175 chtableori<-rbind(sc1,sc2)
177 interori<-MyChiSq(chtableori,sc,TT)/TT#chisq.test(chtableori)$statistic#/TT
178 txt <- paste(txt, ' - interori : ',interori)
179 #pp('interori',interori)
186 txt <- paste(txt, 'N1:', N1,'-N2:',N2)
192 if(cl[compte]==clnb){
193 chtable[1,l]<-chtable[1,l]-1
194 chtable[2,l]<-chtable[2,l]+1
196 chtable[1,l]<-chtable[1,l]+1
197 chtable[2,l]<-chtable[2,l]-1
199 interswitch<-MyChiSq(chtable,sc,TT)/TT#chisq.test(chtable)$statistic/TT
200 ws<-interori-interswitch
203 interori<-interswitch
204 if(cl[compte]==clnb){
208 listsub[[it]]<-append(listsub[[it]],compte)
213 listsub[[it]]<-append(listsub[[it]],compte)
215 vdelta<-append(vdelta,compte)
220 # for (val in vdelta) {
221 # if (cl[val]==clnb) {
223 # listsub[[it]]<-append(listsub[[it]],val)
226 # listsub[[it]]<-append(listsub[[it]],val)
229 print('###################################')
230 print('longueur < 0')
231 malcl<-length(vdelta)
233 if ((it>1)&&(!is.logical(listsub[[it]]))&&(!is.logical(listsub[[it-1]]))){
234 if (all(listsub[[it]]==listsub[[(it-1)]])){
239 print('###################################')
242 #dtable<-cbind(dtable,'cl'=as.vector(cl))
243 #dtable[,'cl'] <-as.vector(cl)
244 #######################################################################
245 # Fin reclassement des individus #
246 #######################################################################
247 # if (!(length(cl[cl==clnb])==1 || length(cl[cl==clnb+1])==1)) {
248 #t1<-dtable[dtable[,'cl']==clnb,][,-ncol(dtable)]
249 #t2<-dtable[dtable[,'cl']==clnb+1,][,-ncol(dtable)]
250 t1<-dtable[which(cl[,1]==clnb),]#[,-ncol(dtable)]
251 t2<-dtable[which(cl[,1]==clnb+1),]#[,-ncol(dtable)]
252 if (inherits(t1, "numeric")) {
259 if (inherits(t2, "numeric")) {
266 chtable<-rbind(sc1,sc2)
267 inter<-chisq.test(chtable)$statistic/TT
268 pp('last inter',inter)
269 print('=====================')
270 #calcul de la specificite des colonnes
271 mint<-min(nrowt1,nrowt2)
272 maxt<-max(nrowt1,nrowt2)
273 seuil<-round((1.9*(maxt/mint))+1.9,digit=6)
274 #sink('/home/pierre/workspace/iramuteq/dev/findchi2.txt')
275 # print('ATTENTION SEUIL 3,84')
293 #another try#########################################
294 one <- cbind(sc1,sc2)
295 cols <- c(length(which(cl==clnb)), length(which(cl==clnb+1)))
297 colss <- matrix(rep(cols,ncol(dtable)), ncol=2, byrow=TRUE)
299 rows <- cbind(rowSums(zero), rowSums(one))
301 for (m in 1:nrow(rows)) {
302 obs <- t(matrix(c(zero[m,],one[m,]),2,2))
303 E <- outer(rows[m,],cols,'*')/n
304 if ((min(obs[2,])==0) & (min(obs[1,])!=0)) {
306 } else if ((min(obs[1,])==0) & (min(obs[2,])!=0)) {
308 } else if (any(obs < 10)) {
309 chi <- sum((abs(obs - E) - 0.5)^2 / E)
311 chi <- sum(((obs - E)^2)/E)
317 if (obs[2,1] < E[2,1]) {
318 listcol[[clnb]]<-append(listcol[[clnb]],cn[m])
320 } else if (obs[2,2] < E[2,2]) {
321 listcol[[clnb+1]]<-append(listcol[[clnb+1]],cn[m])
326 ######################################################
327 print('resultats elim item')
328 pp(clnb+1,length(listcol[[clnb+1]]))
329 pp(clnb,length(listcol[[clnb]]))
331 pp('ncclnbp',ncclnbp)
332 listrownamedtable<-rownames(dtable)
333 listrownamedtable<-as.integer(listrownamedtable)
334 newcol<-vector(length=nrow(dataori))
335 #remplissage de la nouvelle colonne avec les nouvelles classes
338 newcol[listrownamedtable] <- cl[,1]
339 #recuperation de la classe precedante pour les cases vides
340 print('recuperation classes precedentes')
342 newcol[which(newcol==0)] <- dout[,ncol(dout)][which(newcol==0)]
344 if(!is.null(rowelim)) {
347 tailleclasse<-as.matrix(summary(as.factor(as.character(newcol))))
348 print('tailleclasse')
350 tailleclasse<-as.matrix(tailleclasse[!(rownames(tailleclasse)==0),])
351 plusgrand<-which.max(tailleclasse)
352 #???????????????????????????????????
353 #Si 2 classes ont des effectifs egaux, on prend la premiere de la liste...
354 if (length(plusgrand)>1) {
355 plusgrand<-plusgrand[1]
357 #????????????????????????????????????
359 #constuction du prochain tableau a analyser
360 print('construction tableau suivant')
361 dout<-cbind(dout,newcol)
362 classe<-as.integer(rownames(tailleclasse)[plusgrand])
363 dtable<-dataori[which(newcol==classe),]
364 row.names(dtable)<-rownames(dataori)[which(newcol==classe)]
365 colnames(dtable) <- 1:ncol(dtable)
367 listcolelim<-listcol[[as.integer(classe)]]
368 mother<-listmere[[as.integer(classe)]]
370 listcolelim<-append(listcolelim,listcol[[mother]])
371 mother<-listmere[[mother]]
374 listcolelim<-sort(unique(listcolelim))
375 pp('avant',ncol(dtable))
376 if (!is.logical(listcolelim)){
377 print('elimination colonne')
378 #dtable<-dtable[,-listcolelim]
379 dtable<-dtable[,!(colnames(dtable) %in% listcolelim)]
381 pp('apres',ncol(dtable))
382 #elimination des colonnes ne contenant que des 0
383 print('vire colonne inf 3 dans boucle')
386 dtable<-dtable[,-which(sdt<=3)]
388 #elimination des lignes ne contenant que des 0
389 print('vire ligne vide dans boucle')
390 if (ncol(dtable)==1) {
396 rowelim<-as.integer(rownames(dtable)[which(sdt==0)])
397 print('&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&')
399 print('&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&')
400 dtable<-dtable[-which(sdt==0),]
405 res <- list(n1 = dout, list_mere = listmere, list_fille = list_fille)