diff --git a/meldata.cc b/meldata.cc
index 98f3cac5ccc69080ee9a57ca40e40c9916c3bc1e..9af5137c077af77d04303c3cd490f512e83af09a 100644
--- a/meldata.cc
+++ b/meldata.cc
@@ -155,6 +155,8 @@ namespace Melodic{
{
setup_misc();
numfiles = (int)opts.inputfname.value().size();
+ if((numfiles > 1) && (opts.approach.value()==string("defl") || opts.approach.value()==string("symm")))
+ opts.approach.set_T("tica");
Matrix alldat, tmpData;
alldat = process_file(opts.inputfname.value().at(0), numfiles) / numfiles;
@@ -324,6 +326,10 @@ namespace Melodic{
save4D(stdNoisei,string("Noise_stddev_inv"));
}
+ //Output T- & S-modes
+ save_Tmodes();
+ save_Smodes();
+
//Output mixMatrix
if(mixMatrix.Storage()>0){
saveascii(mixMatrix, opts.outputfname.value() + "_mix");
diff --git a/meldata.h b/meldata.h
index 4452240d9d8126ed2c6b58db0c72198dbbdcd4b3..109cd975fcd702a65fb35b40049b756709132a9f 100644
--- a/meldata.h
+++ b/meldata.h
@@ -77,10 +77,22 @@ namespace Melodic{
inline vector<Matrix>& get_Smodes() {return Smodes;}
inline Matrix& get_Smodes(int what) {return Smodes.at(what);}
inline void add_Smodes(Matrix& Arg) {Smodes.push_back(Arg);}
+ inline void save_Smodes(){
+ Matrix tmp = Smodes.at(0);
+ for(unsigned int ctr = 1; ctr < Smodes.size(); ctr++)
+ tmp |= Smodes.at(ctr);
+ saveascii(tmp,opts.outputfname.value() + "_Smodes");
+ }
inline vector<Matrix>& get_Tmodes() {return Tmodes;}
inline Matrix& get_Tmodes(int what) {return Tmodes.at(what);}
- inline void add_Tmodes(Matrix& Arg) {Tmodes.push_back(Arg);}
+ inline void add_Tmodes(Matrix& Arg) {Tmodes.push_back(Arg);}
+ inline void save_Tmodes(){
+ Matrix tmp = Tmodes.at(0);
+ for(unsigned int ctr = 1; ctr < Tmodes.size(); ctr++)
+ tmp |= Tmodes.at(ctr);
+ saveascii(tmp,opts.outputfname.value() + "_Tmodes");
+ }
void set_TSmode();
diff --git a/melica.cc b/melica.cc
index d79b1e2bf5f721aabc0eb876aeed7dd763ba6157..8d7b3351489b583bd743ab66fd66d4152fcfcc3a 100644
--- a/melica.cc
+++ b/melica.cc
@@ -90,8 +90,9 @@ namespace Melodic{
// rank 1 approximation in the case of multiple input files
- if(melodat.get_numfiles() > 1)
+ if(opts.approach.value() == string("tica"))
{
+ message(" Rank-1 approximation of the time courses; ");
Matrix temp(melodat.get_dewhite() * redUMM);
temp = melodat.expand_dimred(temp);
temp = krapprox(temp,int(temp.Nrows()/melodat.get_numfiles()));
@@ -460,7 +461,10 @@ write_ascii_matrix("Xim",Data.Row(1));
//switch to the chosen method
// cout << endl << "Dim = " << dim << endl << "Samples = " << samples << endl;
- if(opts.approach.value()==string("symm")){
+ if(opts.approach.value()==string("symm") ||
+ opts.approach.value()==string("tica") ||
+ opts.approach.value()==string("parafac") ||
+ opts.approach.value()==string("concat")){
ica_fastica_symm(Data);}
if(opts.approach.value()==string("defl")){
ica_fastica_defl(Data);}
@@ -468,6 +472,7 @@ write_ascii_matrix("Xim",Data.Row(1));
ica_jade(Data);}
if(opts.approach.value()==string("maxent")){
ica_maxent(Data);}
+
if(!no_convergence){//calculate the IC
Matrix temp(melodat.get_unmix()*melodat.get_Data());
diff --git a/melodic.cc b/melodic.cc
index 638423ac201d44feaef5029aac26b4909195aaed..3ab18bb6ac25884c3a3d597f459690d7c81ebccd 100644
--- a/melodic.cc
+++ b/melodic.cc
@@ -313,11 +313,13 @@ Matrix mmall(Log& logger, MelodicOptions& opts,
wherelog,ctr,
melodat.get_mask(),
melodat.get_mean(),3);
+ message(" calculating mixture-model fit "<<endl);
mixmod.fit("GGM");
if(opts.output_MMstats.value()){
+ message(" saving probability map:");
melodat.save4D(mixmod.get_probmap(),
- string("stats/probmap_")+num2str(ctr));
+ string("stats/probmap_")+num2str(ctr));
}
//re-scale spatial maps to mean 0 for nht
@@ -380,6 +382,7 @@ Matrix mmall(Log& logger, MelodicOptions& opts,
<< " Means : " << mixmod.get_means() << endl
<< " Vars. : " << mixmod.get_vars() << endl
<< " Prop. : " << mixmod.get_pi() << endl << endl;
+ message(" saving thresholded Z-stats image:");
melodat.save4D(mixmod.get_threshmaps(),
string("stats/thresh_zstat")+num2str(ctr));
}
diff --git a/melodic.h b/melodic.h
index 868e8cfb1046b8d355240f7ae3bb6b8180729332..88bb55c6dc12561501ee52f7e7b8699014702afe 100644
--- a/melodic.h
+++ b/melodic.h
@@ -32,7 +32,7 @@
namespace Melodic{
-const string version = "3.02";
+const string version = "3.0 alpha 1";
}
diff --git a/meloptions.cc b/meloptions.cc
index 5be925e4a87b16c0dbb05241994aa67dfcfc9655..f6633d75e45219fafab43c79c0b4d5ddc8133afb 100644
--- a/meloptions.cc
+++ b/meloptions.cc
@@ -80,7 +80,7 @@ MelodicOptions* MelodicOptions::gopt = NULL;
}
if (approach.value() != "symm" && approach.value() != "defl" &&
approach.value() != "jade" && approach.value() != "maxent" &&
- approach.value() != "tica" && approach.value() != "parafac"){
+ approach.value() != "tica" && approach.value() != "concat"){
cerr << "ERROR:: unknown approach \n\n";
print_usage(argc,argv);
exit(2);
diff --git a/meloptions.h b/meloptions.h
index 85ceb1db12a370d786b3b8c98edb5f6d79d389f0..fc8bcfff4bc7a505c54b79a543d2ce7c2479d88c 100644
--- a/meloptions.h
+++ b/meloptions.h
@@ -171,7 +171,7 @@ class MelodicOptions {
string("numer of IC's to extract (for deflation approach)"),
false, requires_argument),
approach(string("-a,--approach"), string("symm"),
- string("approach for decomposition: defl, symm"),
+ string("approach for decomposition, 2D: defl, symm (default), 3D: tica (default), concat"),
false, requires_argument),
nonlinearity(string("--nl"), string("pow3"),
string("\tnonlinearity: gauss, tanh, pow3, pow4"),
diff --git a/melreport.cc b/melreport.cc
index 47e84b7ce4f2c08243d7c9e96b45b08432985255..6890041fbba0557242799cc9e1c341a6d4893845 100644
--- a/melreport.cc
+++ b/melreport.cc
@@ -146,10 +146,10 @@ namespace Melodic{
{//plot subject mode
Matrix smode;
- smode = melodat.get_Smodes(0);
+ smode = melodat.get_Smodes(cnum-1);
if(smode.Nrows() > 1){
miscplot newplot;
- newplot.timeseries(smode,
+ newplot.timeseries(smode.t(),
report.appendDir(string("s")+num2str(cnum)+".png"),
string("Subject/Session mode"));
write_ascii_matrix(report.appendDir(string("s")
@@ -161,7 +161,6 @@ namespace Melodic{
+string("s")+num2str(cnum)+".png\"></A><p>" << endl;
}
}//subject mode plot
- cerr << "after subject plot " << endl;
{//plot frequency
miscplot newplot;
RowVector empty(1);