1 #include "seaview.h"
2 #include <ctype.h>
3 #include <FL/Fl_Multiline_Output.H>
4 
5 int nbaagroup=5;
6 int aagroup[30];
7 char aalist[30] = "ACDEFGHIKLMNPQRSTVWY";
8 int TRANSVONLY=0;
9 char infotext[10000];
10 
11 extern void *check_alloc(int nbrelt, int sizelt);
12 
bydefault_aagroups()13 void bydefault_aagroups(){
14   int i;
15 
16   sprintf(aalist, "EDQNHRKILMVAPSGTFYWC");
17   for(i=0;i<7;i++) aagroup[i]=0;
18   for(;i<11;i++) aagroup[i]=1;
19   for(;i<16;i++) aagroup[i]=2;
20   for(;i<18;i++) aagroup[i]=3;
21   for(;i<20;i++) aagroup[i]=4;
22 }
23 
24 
25 
26 /* basecomp */
27 /* Write in comp the composition (percent of 4 bases) of sequence seq */
28 
basecomp(char * seq,double * comp)29 int basecomp(char* seq, double* comp){
30   int i, lg=0;
31 
32   for(i=0;i<4;i++) comp[i]=0.;
33   for(i=0;seq[i];i++){
34     if(seq[i]=='A') comp[0]++;
35     else if(seq[i]=='C') comp[1]++;
36     else if(seq[i]=='G') comp[2]++;
37     else if(seq[i]=='T') comp[3]++;
38     else lg--;
39     lg++;
40   }
41   if(lg<=0) return 0;
42   for(i=0;i<4;i++) comp[i]=comp[i]*100/lg;
43   return lg;
44 }
45 
46 
47 /* aaclass */
48 /* Return the class of aminoacid aa */
49 
aaclass(char aa)50 int aaclass(char aa){
51 
52   int i=0;
53 
54   while(aalist[i]!=aa && aagroup[i]!=-1) i++;
55   if(aagroup[i]==-1) return 5;
56   return(aagroup[i]);
57 }
58 
59 
60 /* aacomp */
61 /* Compute aminoacid composition of proteic sequence seq. */
62 /* The whole composition (percent of 20 aminoacid) is written in comp, */
63 /* ordered as in aalist. */
64 /* Percent of each aminoacid group are written in compgroup. */
65 
aacomp(char * seq,double * comp,double * compgroup)66 int aacomp(char* seq, double* comp, double* compgroup){
67   int i, j, lg=0, nbaa;
68   int c;
69 
70   nbaa=(int)strlen(aalist);
71   for(i=0;i<nbaa;i++) comp[i]=0.;
72   for(i=0;i<nbaagroup;i++) compgroup[i]=0.;
73   for(i=0;seq[i];i++){
74     for(j=0;j<nbaa;j++){
75       if(seq[i]==aalist[j]) { comp[j]++; break; }
76     }
77     if(j==nbaa) {
78       lg--;
79       c = -1;
80       }
81     else {
82       c=aaclass(seq[i]);
83       }
84     if(c!=-1) compgroup[c]++;
85     lg++;
86   }
87   if(lg<=0) return 0;
88   for(i=0;i<nbaa;i++) comp[i]=comp[i]*100/lg;
89   for(i=0;i<nbaagroup;i++) compgroup[i]=compgroup[i]*100/lg;
90   return lg;
91 }
92 
93 
94 
95 /* variable */
96 /* Returns the number of variable sites in sequences sequ */
97 
variable(char ** sequ,int nb)98 int variable(char** sequ, int nb){
99   int i, j, var;
100 
101   var=strlen(sequ[0]);
102   for(i=0;sequ[0][i];i++){
103     for(j=1;j<nb;j++){
104       if(sequ[j][i]!=sequ[0][i]) break;
105     }
106     if(j==nb) var--;
107   }
108   return var;
109 }
110 
111 
112 
113 /* informative */
114 /* Returns the number of informative sites in sequences sequ */
115 
informative(char ** sequ,int nb)116 static int informative(char** sequ, int nb){
117 
118   int i, j, k, info=0, siteok=0;
119   char firstpair;
120 
121 
122   for(i=0;sequ[0][i];i++){
123     firstpair=0;
124     siteok=0;
125     for(j=0;j<nb-1;j++){
126       for(k=j+1;k<nb;k++){
127 	if(sequ[j][i]==sequ[k][i]){
128  	  if(firstpair){
129 	    if(sequ[j][i]!=firstpair) { info++; siteok=1; break; }
130 	  }
131 	  else
132 	    firstpair=sequ[j][i];
133         }
134       }
135       if (siteok) break;
136     }
137   }
138 
139   return info;
140 }
141 
142 
143 /* refresh */
144 /* Remove all sites with gaps only from sequence alignment seq. */
145 /* If option==0, all sites with at least one gap are removed. */
146 
refresh(char ** seq,int nbseq,int option,int prot)147 void refresh(char** seq, int nbseq, int option, int prot)
148 {
149 
150   int lgseq, l=-1, drapeau, i, j, k;
151   char **seqref ;
152 
153   lgseq=(int)strlen(seq[0]);
154   seqref = (char**)check_alloc(nbseq, sizeof(char*));
155   for(i=0;i<nbseq;i++)
156     seqref[i]=(char*)check_alloc(lgseq+1, sizeof(char));
157 
158   if (option==0){ /* no site with at least one gap */
159     for(i=0;i<lgseq;i++){
160       drapeau=0;
161       for(j=0;j<nbseq;j++){
162 	if (!prot  && seq[j][i]!='A' && seq[j][i]!='C' && seq[j][i]!='G' && seq[j][i]!='T' && (!TRANSVONLY || (seq[j][i]!='R' && seq[j][i]!='Y'))) drapeau=1;
163 	if (prot && (seq[j][i]=='X' || seq[j][i]=='-')) drapeau=1;
164       }
165       if (drapeau==0){
166 	l++;
167 	for(k=0;k<nbseq;k++) *(seqref[k]+l)=*(seq[k]+i);
168       }
169     }
170     for(i=0;i<nbseq;i++)
171       for (j=l+1;j<lgseq;j++)
172         *(seqref[i]+j)='\0';
173   }
174   else{ /* no site with only gaps */
175     for(i=0;i<lgseq;i++){
176       drapeau=0;
177       for(j=0;j<nbseq;j++){
178 	if (!prot && ((*(seq[j]+i)=='A' || *(seq[j]+i)=='C' || *(seq[j]+i)=='T' || *(seq[j]+i)=='G') || (TRANSVONLY && (*(seq[j]+i)=='R' || *(seq[j]+i)=='Y')))) {
179 	  drapeau=1;
180 	  break;
181 	}
182 	if (prot && (*(seq[j]+i)!='X' || *(seq[j]+i)!='-')) {
183 	  drapeau=1;
184 	  break;
185 	}
186       }
187       if (drapeau==1){
188 	l++;
189 	for(k=0;k<nbseq;k++){
190 	  *(seqref[k]+l)=*(seq[k]+i);
191 	  *(seqref[k]+l+1)=0;
192         }
193       }
194     }
195   }
196   for (i=0;i<nbseq;i++)
197     for (j=0;j<lgseq;j++)
198       *(seq[i]+j)=*(seqref[i]+j);
199 
200   for(i=0;i<nbseq;i++) free(seqref[i]);
201   free(seqref);
202 }
203 
204 
reg_and_pos(SEA_VIEW * view,int * selected,list_segments * active)205 void reg_and_pos(SEA_VIEW *view, int *selected, list_segments* active)
206 {
207 
208   int sel1=0, sel2=0, sel3=0, i /*, cur_reg_deb=0*/;
209   int lgseq = view->seq_length;
210   list_segments start_end = {1, lgseq, NULL};
211   if (active == NULL) active = &start_end;
212 
213   /*if(view->protein){
214     sel1=sel2=sel3=1;
215   }
216   else{
217     if(GetStatus(po1i)) sel1=1;
218     if(GetStatus(po2i)) sel2=1;
219     if(GetStatus(po3i)) sel3=1;
220   }*/
221   sel1=sel2=sel3=1;
222 
223   if(sel1 && sel2 && sel3){
224     while(active) {
225       for(i = active->debut - 1; i < active->fin; i++) selected[i]=1;
226       active = active->next;
227       }
228   }
229   else{
230     /*if(sel1 && sel2 && !sel3){
231       for(i=0;i<lgseq;i++){
232         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
233 	j=i-cur_reg_deb;
234         if(regions[i] && (j%3==0 || j%3==1)) selected[i]=1; else selected[i]=0;
235       }
236     }
237     if(sel1 && !sel2 && sel3){
238       for(i=0;i<lgseq;i++){
239         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
240 	j=i-cur_reg_deb;
241 	if(regions[i] && (j%3==0 || j%3==2)) selected[i]=1; else selected[i]=0;
242       }
243     }
244     if(!sel1 && sel2 && sel3){
245       for(i=0;i<lgseq;i++){
246         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
247 	j=i-cur_reg_deb;
248 	if(regions[i] && (j%3==1 || j%3==2)) selected[i]=1; else selected[i]=0;
249       }
250     }
251     if(sel1 && !sel2 && !sel3){
252       for(i=0;i<lgseq;i++){
253         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
254 	j=i-cur_reg_deb;
255 	if(regions[i] && (j%3==0)) selected[i]=1; else selected[i]=0;
256       }
257     }
258     if(!sel1 && sel2 && !sel3){
259       for(i=0;i<lgseq;i++){
260         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
261 	j=i-cur_reg_deb;
262 	if(regions[i] && (j%3==1)) selected[i]=1; else selected[i]=0;
263       }
264     }
265     if(!sel1 && !sel2 && sel3){
266       for(i=0;i<lgseq;i++){
267         if(regions[i] && !regions[i-1]) cur_reg_deb=i;
268 	j=i-cur_reg_deb;
269  	if(regions[i] && (j%3==2)) selected[i]=1; else selected[i]=0;
270       }
271     }*/
272     if(!sel1 && !sel2 && !sel3){
273       for(i=0;i<lgseq;i++) selected[i]=0;
274     }
275   }
276   return;
277 }
278 
279 
280 /* freq_obs */
281 /* Write at address freq observed frequencies of 16 di-nucleotides XY */
282 /* (X= A,C,G,T  ,  Y=A,C,G,T) X and Y being homologous nucleotides of sequences */
283 /* seq1 and seq2. Alphabetic order is used : freq[0]=AA frequency, freq[1]=AC, */
284 /* ..., freq[15]=TT. */
freq_obs(char * seq1,char * seq2,double * freq)285 int freq_obs(char* seq1, char* seq2, double* freq){
286 
287   int i, lgseq, lgseqvrai;
288 
289   if ((int)strlen(seq1)!=(int)strlen(seq2)){
290     printf ("Longueurs inegales.\n");
291     return -1;
292   }
293   for(i=0;i<16;i++) freq[i]=0;
294   lgseq=(int)strlen(seq1);
295   lgseqvrai=lgseq;
296   for(i=0;i<lgseq;i++){
297     switch(seq1[i]){
298       case 'A':
299 	switch(seq2[i]){
300 	  case 'A' : freq[0]++; break;
301 	  case 'C' : freq[1]++; break;
302 	  case 'G' : freq[2]++; break;
303 	  case 'T' : freq[3]++; break;
304 	  default : lgseqvrai --; break;
305 	}
306 	break;
307       case 'C':
308 	switch(seq2[i]){
309 	  case 'A' : freq[4]++; break;
310 	  case 'C' : freq[5]++; break;
311 	  case 'G' : freq[6]++; break;
312 	  case 'T' : freq[7]++; break;
313 	  default : lgseqvrai --; break;
314 	}
315 	break;
316       case 'G':
317 	switch(seq2[i]){
318 	  case 'A' : freq[8]++; break;
319 	  case 'C' : freq[9]++; break;
320 	  case 'G' : freq[10]++; break;
321 	  case 'T' : freq[11]++; break;
322 	  default : lgseqvrai --; break;
323 	}
324 	break;
325       case 'T':
326 	switch(seq2[i]){
327 	  case 'A' : freq[12]++; break;
328 	  case 'C' : freq[13]++; break;
329 	  case 'G' : freq[14]++; break;
330 	  case 'T' : freq[15]++; break;
331 	  default : lgseqvrai --; break;
332 	}
333 	break;
334       default :
335 	lgseqvrai --;
336     }
337   }
338   if(lgseqvrai!=0){
339     for(i=0;i<16;i++) freq[i]/=lgseqvrai;
340     return 1;
341   }
342   else return 0;
343 }
344 
345 
save_stat_cb(Fl_Widget * ob,void * data)346 static void save_stat_cb(Fl_Widget *ob, void *data)
347 {
348   Fl_Multiline_Output *ml = (Fl_Multiline_Output*)data;
349   Fl_Native_File_Chooser fnfc;
350   fnfc.title("Pick a text file");
351   fnfc.type(Fl_Native_File_Chooser::BROWSE_SAVE_FILE);
352   fnfc.filter("Text\t*.txt\n");
353   fnfc.preset_file("statistics.txt");
354   fnfc.options(Fl_Native_File_Chooser::SAVEAS_CONFIRM | fnfc.options());
355   if ( fnfc.show() != 0 ) return;
356   FILE *out = fopen(fnfc.filename(), "w");
357   fprintf(out, "Alignment: %s\n\n", (char*)ml->user_data());
358   fputs(ml->value(), out);
359   fclose(out);
360 }
361 
geni_act(SEA_VIEW * view)362 void geni_act(SEA_VIEW *view)
363 {
364   int i, j, k, l, ii, nb, lgtot, lg, lgvar, lginfo, *lgcomp, totlgcomp, nc=0;
365   int *selected, nbaa, nbcp;
366   double **comp, *compall, compcomp[4], **compgroup, *compgroupall, min[4], max[4];
367   double P, Psum=0., Q, Qsum=0., Pmin=0, Pmax=0, Qmin=0, Qmax=0, titvratio, titvmin, titvmax, freq[16], mat_obs[16], *ti, *tv;
368   char **genseq, **genname, *text, *minname[4], *maxname[4], nt[4]={'A', 'C', 'G', 'T'};
369   const char *titvminname1 = "", *titvminname2 = "";
370   const char *titvmaxname1 = "", *titvmaxname2 = "";
371   int lgseq = view->seq_length;
372   int nbseq = view->tot_seqs;
373 
374   bydefault_aagroups();
375 
376   text=(char*)check_alloc(10000, sizeof(char));
377 
378   nb=lgtot=0;
379   selected=(int*)check_alloc(lgseq, sizeof(int));
380   reg_and_pos(view, selected, view->active_region ? view->active_region->list : NULL);
381   nbaa=(int)strlen(aalist);
382 
383   /* SET UP DATA */
384 
385   for(i=0;i<lgseq;i++) if (selected[i]) lgtot++;
386   if(lgtot==0){
387     fl_alert("No selected sites");
388     free(selected); free(text);
389     return;
390   }
391   for(i=0;i<nbseq;i++) if (view->tot_sel_seqs == 0 || view->sel_seqs[i]) nb++;
392   genseq=(char**)check_alloc(nb, sizeof(char*));
393   genname=(char**)check_alloc(nb, sizeof(char*));
394   for(i=0;i<nb;i++)
395     genseq[i]=(char*)check_alloc(lgtot+1, sizeof(char));
396   j=0;
397   for(i=0;i<nbseq;i++){
398     if(view->tot_sel_seqs == 0 || view->sel_seqs[i]){
399       genname[j]=view->seqname[i];
400       l=0;
401       for(k=0;k<lgseq;k++) {
402 	if(selected[k]) genseq[j][l++] = toupper(view->sequence[i][k]);
403 	}
404       j++;
405     }
406   }
407 
408 
409   ti=(double*)check_alloc(nb*(nb-1)/2, sizeof(double));
410   tv=(double*)check_alloc(nb*(nb-1)/2, sizeof(double));
411 
412 
413   /* ALL SITES : COMPOSITION */
414 
415   lgcomp=(int*)check_alloc(nb, sizeof(int));
416   comp=(double**)check_alloc(nb, sizeof(double*));
417   for(i=0;i<nb;i++)
418     comp[i]=(double*)check_alloc(view->protein?nbaa:4, sizeof(double));
419   compall=(double*)check_alloc(view->protein?nbaa:4, sizeof(double));
420 
421   if(view->protein){
422     compgroup=(double**)check_alloc(nb, sizeof(double*));
423     for(i=0;i<nb;i++)
424       compgroup[i]=(double*)check_alloc(nbaagroup, sizeof(double));
425     compgroupall=(double*)check_alloc(nbaagroup, sizeof(double));
426 
427     totlgcomp=0;
428     for(i=0;i<nb;i++){
429       lgcomp[i]=aacomp(genseq[i], comp[i], compgroup[i]);
430       totlgcomp+=lgcomp[i];
431     }
432     for(i=0;i<nbaa;i++){
433       compall[i]=0.;
434       for(j=0;j<nb;j++)
435 	compall[i]+=comp[j][i]*lgcomp[j];
436       compall[i]/=totlgcomp;
437     }
438     for(i=0;i<nbaagroup;i++){
439       compgroupall[i]=0.;
440       for(j=0;j<nb;j++)
441 	compgroupall[i]+=compgroup[j][i]*lgcomp[j];
442       compgroupall[i]/=totlgcomp;
443     }
444   }
445   else{
446     totlgcomp=0;
447     for(i=0;i<nb;i++){
448       lgcomp[i]=basecomp(genseq[i], comp[i]);
449       totlgcomp+=lgcomp[i];
450     }
451     for(i=0;i<4;i++){
452       compall[i]=0.;
453       for(j=0;j<nb;j++)
454 	compall[i]+=comp[j][i]*lgcomp[j];
455       compall[i]/=totlgcomp;
456 
457     }
458   }
459 
460 
461   /* COMPLETE SITES ONLY */
462 
463   refresh(genseq, nb, 0, view->protein);
464   lg=(int)strlen(genseq[0]);
465   lgvar=variable(genseq, nb);
466   lginfo=informative(genseq, nb);
467 
468   /* COMPOSITION */
469 
470   if(!view->protein && lg!=0){
471     for(i=0;i<nb;i++)
472       lgcomp[i]=basecomp(genseq[i], comp[i]);
473     for(i=0;i<4;i++){
474       compcomp[i]=0;
475       for(j=0;j<nb;j++) compcomp[i]+=comp[j][i];
476       compcomp[i]/=nb;
477       min[i]=max[i]=comp[0][i]; minname[i]=maxname[i]=genname[0];
478       for(j=0;j<nb;j++){
479 	if(comp[j][i]<min[i]){ min[i]=comp[j][i]; minname[i]=genname[j]; }
480 	if(comp[j][i]>max[i]){ max[i]=comp[j][i]; maxname[i]=genname[j]; }
481       }
482     }
483   }
484 
485   /* CHANGES */
486 
487   if(!view->protein && lg!=0){
488 
489     for(ii=0;ii<16;ii++) mat_obs[ii]=0.;
490     k=0; nbcp=0;
491     titvmin=100000.; titvmax=0.;
492     for(i=0;i<nb-1;i++){
493       for(j=i+1;j<nb;j++){
494         freq_obs(genseq[i], genseq[j], freq);
495         P=(freq[2]+freq[7]+freq[8]+freq[13])*lg;
496         Q=(freq[1]+freq[3]+freq[4]+freq[6]+freq[9]+freq[11]+freq[12]+freq[14])*lg;
497         for(ii=0;ii<16;ii++) mat_obs[ii]+=freq[ii];
498 
499 	ti[nc]=P;
500 	tv[nc]=Q;
501 	nc++;
502 
503 
504 	Psum+=P;
505 	Qsum+=Q;
506         if(Q>0.1){
507 	  titvratio=P/Q;
508           if(titvratio<titvmin){
509 	    titvmin=titvratio;
510 	    titvminname1=genname[i];
511 	    titvminname2=genname[j];
512 	    Pmin=P; Qmin=Q;
513 	  }
514       	  if(titvratio>titvmax){
515 	    titvmax=titvratio;
516 	    titvmaxname1=genname[i];
517 	    titvmaxname2=genname[j];
518 	    Pmax=P; Qmax=Q;
519 	  }
520 	  nbcp++;
521         }
522         k++;
523       }
524     }
525 
526     for(ii=0;ii<16;ii++) mat_obs[ii]/=(nb*(nb-1)/(2*100.));
527 
528     /*
529      printf("observed matrix :\n");
530      printf("AA : %.1f  CC : %.1f  GG : %.1f  TT : %.1f\n", mat_obs[0], mat_obs[5], mat_obs[10], mat_obs[15]);
531      printf("AG : %.1f  CT : %.1f\n", mat_obs[2]+mat_obs[8], mat_obs[7]+mat_obs[13]);
532      printf("AC : %.1f  AT : %.1f  CG : %.1f  GT : %.1f\n", mat_obs[1]+mat_obs[4], mat_obs[3]+mat_obs[12], mat_obs[6]+mat_obs[9], mat_obs[11]+mat_obs[14]);
533      */
534 
535     /*
536      titvf=fopen("TITV", "w");
537      fprintf(titvf, "Ti\tTv\n");
538      for(i=0;i<nc;i++) fprintf(titvf, "%f\t%f\n", ti[i], tv[i]);
539      fclose(titvf);
540      */
541 
542   }
543 
544 
545   /* PAIRWISE SUBSTITUTION MATRIX */
546 
547   if(!view->protein && nb==2 && lg!=0)
548     freq_obs(genseq[0], genseq[1], freq);
549 
550 
551   /* DISPLAY INTO NEW WINDOW */
552   char *p;
553   Fl_Window *win = new Fl_Window(600, 500, "Statistics");
554   Fl_Multiline_Output *ml = new Fl_Multiline_Output(2, 2, win->w()-4, win->h()-24);
555   ml->textfont(FL_COURIER);
556   ml->user_data(view->masename);
557   win->resizable(ml);
558   Fl_Button *b = new Fl_Button(win->w()/2 - 25 , win->h() - 23, 50, 20, "Save");
559   b->callback(save_stat_cb, ml);
560   Fl_Return_Button *ok = new Fl_Return_Button(win->w() - 55 , win->h() - 23, 50, 20, "OK");
561 
562   if(view->tot_sel_seqs == 0) sprintf(text, "%d species\n", nb);
563   else {
564     sprintf(text, "%d selected species\n", nb);
565   }
566   if(view->tot_sel_seqs > 0 && nb<=20){
567     for(i=0;i<view->tot_seqs;i++){
568       if(!view->sel_seqs[i]) continue;
569       strcat(text, view->seqname[i]); strcat(text, ", ");
570       if((i+1)%4==0 || i==nb-1) strcat(text, "\n");
571     }
572     strcat(text, "\n");
573   }
574   p = text + strlen(text);
575   sprintf(p, "\n%d selected sites\n", lgtot);
576   p += strlen(p);
577   if(view->protein)
578     sprintf(p, "    including    %d complete (no gaps, no X)\n", lg);
579   else
580     sprintf(p, "    including    %d complete (no gaps, no N)\n", lg);
581   p += strlen(p);
582   if(nb>=2 && lg!=0){
583     sprintf(p, "    including    %d variable (%.1f%% of complete)\n", lgvar, (double)lgvar*100/lg);
584     p += strlen(p);
585     if(nb>=4){
586       sprintf(p, "    including    %d informative (%.1f%% of complete)\n", lginfo, (double)lginfo*100/lg);
587       p += strlen(p);
588     }
589   }
590   if(view->protein){
591     int k;
592     sprintf(p, "\nAMINOACID COMPOSITION (all sites) :\n");
593     p += strlen(p);
594     j=0;
595     for(i=0;i<nbaagroup;i++){
596       sprintf(p, "\n\nGroup %d : %.1f%%", i+1, compgroupall[i]);
597       p += strlen(p);
598       k = 0;
599       while(aagroup[j]==i){
600 	if(k++ % 3 == 0) {
601 	  sprintf(p, "\n    "); p += strlen(p);
602 	  }
603 	sprintf(p, "%c : %.1f     ", aalist[j], compall[j]);
604 	p += strlen(p);
605 	j++;
606       }
607       if(k % 3 != 0) sprintf(p++, "\n");
608     }
609     if(j!=nbaa){
610       sprintf(p, "\n\nOthers : \n");
611       p += strlen(p);
612       while(aalist[j]){
613 	sprintf(p, "%c : %.1f     ", aalist[j], compall[j]);
614 	p += strlen(p);
615 	sprintf(p++, "\n");
616 	j++;
617       }
618     }
619   }
620   else{
621     sprintf(p, "\nBASE COMPOSITION :\n");
622     p += strlen(p);
623     sprintf(p, "     All sites  :\t %.1f%% A   %.1f%% C   %.1f%% G   %.1f%% T\n", compall[0], compall[1], compall[2], compall[3]);
624     p += strlen(p);
625     if(nb>1 && lg!=0){
626       sprintf(p, "     Complete sites only :\n");
627       p += strlen(p);
628       sprintf(p, "                                Minimum                   Maximum:\n");
629       p += strlen(p);
630       for(i=0;i<4;i++){
631         sprintf(p, "     %c : %.1f %%    ", nt[i], compcomp[i]);
632 	p += strlen(p);
633         sprintf(p, "%.1f%% (%s)     ", min[i], minname[i]);
634  	p += strlen(p);
635         sprintf(p, "%.1f%% (%s)\n", max[i], maxname[i]);
636 	p += strlen(p);
637       }
638     }
639   }
640 
641   if(!view->protein && nb>1 && lg!=0){
642     sprintf(p, "\nOBSERVED CHANGES (complete sites)\n"); p += strlen(p);
643     double PQratio = Qsum ? (double)Psum/Qsum : 0;
644     if(nb>2)
645       sprintf(p, "Transition/transversion ratio : %.3f (mean over all sequence pairs)\n", PQratio);
646     else
647       sprintf(p, "Transition/transversion ratio : %.3f\n", PQratio);
648     p += strlen(p);
649     if(nb>2){
650       sprintf(p, "Minimum : %.3f (%s - %s, %d transitions, %d transversions)\n", titvmin, titvminname1, titvminname2, (int)Pmin, (int)Qmin);
651       p += strlen(p);
652       sprintf(p, "Maximum : %.3f (%s - %s, %d transitions, %d transversions)\n", titvmax, titvmaxname1, titvmaxname2, (int)Pmax, (int)Qmax);
653       //p += strlen(p);
654     }
655     else {
656       sprintf(p, "Observed nucleotide pairs :\n"); p += strlen(p);
657       sprintf(p, "%35s%s\n"
658 	         "%36s A     C     G     T\n", "", genname[1], ""); p += strlen(p);
659       for(i=0;i<4;i++){
660 	sprintf(p, "%30s %c", i==1 ? genname[0] : "", nt[i]); p += strlen(p);
661    	for(j=0;j<4;j++){
662 	  sprintf(p, "%6d", (int)(freq[4*i+j]*lgtot)); p += strlen(p);
663 	}
664 	*p++ = '\n';
665       }
666     }
667   }
668 
669   ml->value(text);
670   win->end();
671   win->show();
672   while (TRUE) {
673     Fl_Widget* obj = Fl::readqueue();
674     if (!obj) Fl::wait();
675     else if (obj == ok || obj == win) {
676       delete win;
677       break;
678     }
679   }
680 
681   free(selected); free(text);
682   for(i=0;i<nb;i++) {free(genseq[i]); free(comp[i]);}
683   free(genseq); free(genname); free(comp);
684   free(lgcomp); free(compall);
685   if(view->protein){
686     for(i=0;i<nb;i++) free(compgroup[i]);
687     free(compgroup);
688     free(compgroupall);
689   }
690 }
691 
stats_callback(Fl_Widget * obj,void * data)692 void stats_callback(Fl_Widget *obj, void *data)
693 {
694   geni_act((SEA_VIEW*)obj->user_data());
695 }
696