Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Shape/Scalar Correlation #2168

Merged
merged 24 commits into from
Nov 28, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
24 commits
Select commit Hold shift + click to select a range
d4cd3a4
Add interface and infrastrucutre for shape scalar correlation.
akenmorris Oct 10, 2023
802cf36
Add mbpls==1.0.4
akenmorris Oct 10, 2023
b1a9e0a
Run 2 block pls
akenmorris Oct 19, 2023
b68df4e
Correctly mark job as complete.
akenmorris Oct 30, 2023
9a1e819
Run mbpls, generate MSE plot image and return.
akenmorris Oct 30, 2023
d1cbe3b
Fix eigen blocking for 2 block PLS
akenmorris Oct 31, 2023
7154cc4
Fix 2blockPLS, formulation was incorrect.
akenmorris Nov 3, 2023
b52e7dc
Cleanup
akenmorris Nov 3, 2023
407675f
Adding structure to find number of components for 2bpls
akenmorris Nov 8, 2023
fd5f9e5
Add support for determining the number of components to use.
akenmorris Nov 9, 2023
1920992
Adding UI for shape/scalar PCA stuff
akenmorris Nov 15, 2023
1fca3d1
Some cleanup
akenmorris Nov 17, 2023
3531f6d
Convert all EigenXf to EigenXd.
akenmorris Nov 17, 2023
216ea05
Cleanup ParticleShapeStatistics
akenmorris Nov 18, 2023
6632c1e
More cleanup
akenmorris Nov 18, 2023
deeebdd
Fix stats for position + scalar
akenmorris Nov 18, 2023
410c29c
Extrat shape from combined PCA.
akenmorris Nov 18, 2023
71a99a1
PCA on shape+scalars
akenmorris Nov 20, 2023
f87ebe3
New UI for shape/scalar stuff.
akenmorris Nov 20, 2023
181baf4
Speed up 2bpls prediction.
akenmorris Nov 20, 2023
bed1237
Scalars only PCA modes
akenmorris Nov 21, 2023
7e2cf17
Remaining pieces for pca shape scalar
akenmorris Nov 21, 2023
886f8ef
Revert feature_map_example.xlsx
akenmorris Nov 27, 2023
d841cfa
Fixes for MS compiler
akenmorris Nov 27, 2023
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
56 changes: 28 additions & 28 deletions Libs/Analyze/Analyze.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@ namespace shapeworks {

//---------------------------------------------------------------------------
static json get_eigen_vectors(ParticleShapeStatistics* stats) {
auto values = stats->Eigenvectors();
auto values = stats->get_eigen_vectors();

std::vector<double> vals;
for (size_t i = values.cols() - 1, ii = 0; i > 0; i--, ii++) {
Expand Down Expand Up @@ -173,8 +173,8 @@ void Analyze::run_offline_analysis(std::string outfile, float range, float steps
double increment = range / half_steps;

std::vector<double> eigen_vals;
for (int i = stats_.Eigenvalues().size() - 1; i > 0; i--) {
eigen_vals.push_back(stats_.Eigenvalues()[i]);
for (int i = stats_.get_eigen_values().size() - 1; i > 0; i--) {
eigen_vals.push_back(stats_.get_eigen_values()[i]);
}
double sum = std::accumulate(eigen_vals.begin(), eigen_vals.end(), 0.0);

Expand Down Expand Up @@ -213,7 +213,7 @@ void Analyze::run_offline_analysis(std::string outfile, float range, float steps
// export modes
std::vector<json> modes;
for (int mode = 0; mode < num_modes; mode++) {
unsigned int m = stats_.Eigenvectors().cols() - (mode + 1);
unsigned int m = stats_.get_eigen_vectors().cols() - (mode + 1);
json jmode;
jmode["mode"] = mode + 1;
double eigen_value = eigen_vals[mode];
Expand All @@ -231,8 +231,8 @@ void Analyze::run_offline_analysis(std::string outfile, float range, float steps
jmode["cumulative_explained_variance"] = cumulative_explained_variance;
SW_LOG("explained_variance [{}]: {:.2f}", mode, explained_variance);
SW_LOG("cumulative_explained_variance [{}]: {:.2f}", mode, cumulative_explained_variance);

double lambda = sqrt(stats_.Eigenvalues()[m]);
double lambda = sqrt(stats_.get_eigen_values()[m]);

std::vector<json> jmodes;
for (int i = 1; i <= half_steps; i++) {
Expand Down Expand Up @@ -278,7 +278,7 @@ void Analyze::run_offline_analysis(std::string outfile, float range, float steps
}

j["eigen_vectors"] = get_eigen_vectors(&stats_);
j["eigen_values"] = stats_.Eigenvalues();
j["eigen_values"] = stats_.get_eigen_values();
j["modes"] = modes;
j["charts"] = create_charts(&stats_);

Expand Down Expand Up @@ -322,7 +322,7 @@ Particles Analyze::get_mean_shape_points() {
return Particles();
}

return convert_from_combined(stats_.Mean());
return convert_from_combined(stats_.get_mean());
}

//---------------------------------------------------------------------------
Expand All @@ -336,30 +336,30 @@ ShapeHandle Analyze::get_mean_shape() {

//---------------------------------------------------------------------------
Particles Analyze::get_shape_points(int mode, double value) {
if (!compute_stats() || stats_.Eigenvectors().size() <= 1) {
if (!compute_stats() || stats_.get_eigen_vectors().size() <= 1) {
return Particles();
}
if (mode + 2 > stats_.Eigenvalues().size()) {
mode = stats_.Eigenvalues().size() - 2;
if (mode + 2 > stats_.get_eigen_values().size()) {
mode = stats_.get_eigen_values().size() - 2;
}

unsigned int m = stats_.Eigenvectors().cols() - (mode + 1);

Eigen::VectorXd e = stats_.Eigenvectors().col(m);

double lambda = sqrt(stats_.Eigenvalues()[m]);
unsigned int m = stats_.get_eigen_vectors().cols() - (mode + 1);
Eigen::VectorXd e = stats_.get_eigen_vectors().col(m);
double lambda = sqrt(stats_.get_eigen_values()[m]);

std::vector<double> vals;
for (int i = stats_.Eigenvalues().size() - 1; i > 0; i--) {
vals.push_back(stats_.Eigenvalues()[i]);
for (int i = stats_.get_eigen_values().size() - 1; i > 0; i--) {
vals.push_back(stats_.get_eigen_values()[i]);
}
double sum = std::accumulate(vals.begin(), vals.end(), 0.0);
double cumulation = 0;
for (size_t i = 0; i < mode + 1; ++i) {
cumulation += vals[i];
}

auto temp_shape = stats_.Mean() + (e * (value * lambda));
auto temp_shape = stats_.get_mean() + (e * (value * lambda));

return convert_from_combined(temp_shape);
}
Expand Down Expand Up @@ -489,9 +489,9 @@ bool Analyze::compute_stats() {
return false;
}
}

stats_.ImportPoints(points, group_ids);
stats_.ComputeModes();
stats_.import_points(points, group_ids);
stats_.compute_modes();

stats_ready_ = true;
SW_LOG("Computed stats successfully");
Expand Down Expand Up @@ -522,7 +522,7 @@ Particles Analyze::convert_from_combined(const Eigen::VectorXd& points) {

//---------------------------------------------------------------------------
void Analyze::initialize_mesh_warper() {
int median = stats_.ComputeMedianShape(-32); //-32 = both groups
int median = stats_.compute_median_shape(-32); //-32 = both groups

if (median < 0 || median >= get_shapes().size()) {
SW_ERROR("Unable to set reference mesh, stats returned invalid median index");
Expand Down Expand Up @@ -555,9 +555,9 @@ void Analyze::initialize_mesh_warper() {
int Analyze::get_num_subjects() { return shapes_.size(); }

//---------------------------------------------------------------------------
Eigen::VectorXf Analyze::get_subject_features(int subject, std::string feature_name) {
Eigen::VectorXd Analyze::get_subject_features(int subject, std::string feature_name) {
if (subject >= shapes_.size()) {
return Eigen::VectorXf();
return Eigen::VectorXd();
}

auto shape = shapes_[subject];
Expand Down Expand Up @@ -617,8 +617,8 @@ Particles Analyze::get_group_shape_particles(double ratio) {
if (!compute_stats()) {
return Particles();
}

auto particles = stats_.Group1Mean() + (stats_.GroupDifference() * ratio);
auto particles = stats_.get_group1_mean() + (stats_.get_group_difference() * ratio);

return convert_from_combined(particles);
}
Expand Down
2 changes: 1 addition & 1 deletion Libs/Analyze/Analyze.h
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ class Analyze {

ShapeHandle create_shape_from_points(Particles points);

Eigen::VectorXf get_subject_features(int subject, std::string feature_name);
Eigen::VectorXd get_subject_features(int subject, std::string feature_name);

void set_group_selection(std::string group_name, std::string group1, std::string group2);

Expand Down
12 changes: 6 additions & 6 deletions Libs/Analyze/Shape.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -572,7 +572,7 @@ void Shape::apply_feature_to_points(std::string feature, ImageType::Pointer imag

int num_points = all_locals.size() / 3;

Eigen::VectorXf values(num_points);
Eigen::VectorXd values(num_points);

int idx = 0;
for (int i = 0; i < num_points; ++i) {
Expand Down Expand Up @@ -615,7 +615,7 @@ void Shape::load_feature_from_mesh(std::string feature, MeshHandle mesh) {

int num_points = all_locals.size() / 3;

Eigen::VectorXf values(num_points);
Eigen::VectorXd values(num_points);

vtkDataArray* from_array = from_mesh->GetPointData()->GetArray(feature.c_str());
if (!from_array) {
Expand All @@ -639,10 +639,10 @@ void Shape::load_feature_from_mesh(std::string feature, MeshHandle mesh) {
}

//---------------------------------------------------------------------------
Eigen::VectorXf Shape::get_point_features(std::string feature) {
Eigen::VectorXd Shape::get_point_features(std::string feature) {
auto it = point_features_.find(feature);
if (it == point_features_.end()) {
return Eigen::VectorXf();
return Eigen::VectorXd();
}

return it->second;
Expand Down Expand Up @@ -683,7 +683,7 @@ std::vector<vtkSmartPointer<vtkTransform>> Shape::get_procrustes_transforms() {
}

//---------------------------------------------------------------------------
void Shape::set_point_features(std::string feature, Eigen::VectorXf values) {
void Shape::set_point_features(std::string feature, Eigen::VectorXd values) {
point_features_[feature] = values;

auto group = get_meshes(DisplayMode::Reconstructed);
Expand Down Expand Up @@ -778,7 +778,7 @@ void Shape::load_feature_from_scalar_file(std::string filename, std::string feat
floats.push_back(line);
}

Eigen::VectorXf values(floats.size());
Eigen::VectorXd values(floats.size());
for (int i = 0; i < floats.size(); i++) {
values[i] = floats[i];
}
Expand Down
6 changes: 3 additions & 3 deletions Libs/Analyze/Shape.h
Original file line number Diff line number Diff line change
Expand Up @@ -158,9 +158,9 @@ class Shape {

std::shared_ptr<Image> get_image_volume(std::string image_volume_name);

Eigen::VectorXf get_point_features(std::string feature);
Eigen::VectorXd get_point_features(std::string feature);

void set_point_features(std::string feature, Eigen::VectorXf values);
void set_point_features(std::string feature, Eigen::VectorXd values);

void load_feature_from_scalar_file(std::string filename, std::string feature_name);

Expand Down Expand Up @@ -198,7 +198,7 @@ class Shape {
std::vector<std::string> global_point_filenames_;
std::vector<std::string> local_point_filenames_;

std::map<std::string, Eigen::VectorXf> point_features_;
std::map<std::string, Eigen::VectorXd> point_features_;
Particles particles_;

std::shared_ptr<shapeworks::Subject> subject_;
Expand Down
2 changes: 1 addition & 1 deletion Libs/Analyze/StudioMesh.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ void StudioMesh::apply_feature_map(std::string name, ImageType::Pointer image) {

//---------------------------------------------------------------------------
void StudioMesh::interpolate_scalars_to_mesh(std::string name, Eigen::VectorXd positions,
Eigen::VectorXf scalar_values) {
Eigen::VectorXd scalar_values) {
int num_points = positions.size() / 3;
if (num_points == 0) {
return;
Expand Down
2 changes: 1 addition & 1 deletion Libs/Analyze/StudioMesh.h
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@ class StudioMesh {
void apply_scalars(MeshHandle mesh);

//! Interpolation scalars at positions to this mesh
void interpolate_scalars_to_mesh(std::string name, Eigen::VectorXd positions, Eigen::VectorXf scalar_values);
void interpolate_scalars_to_mesh(std::string name, Eigen::VectorXd positions, Eigen::VectorXd scalar_values);

//! Return the range of largest axis (e.g. 200 for an object that sits in 100x200x100)
double get_largest_dimension_size();
Expand Down
Loading
Loading