sourcesinkdata.h

Go to the documentation of this file.

#ifndef SOURCESINKDATA_H
#define SOURCESINKDATA_H
 
#include <elemental_profile_set.h>
#include "vector"
#include "parameter.h"
#include "observation.h"
#include "contribution.h"
#include "results.h"
#include "progressreporter.h"
#include "fstream"
#include "qjsonobject.h"
#include "cmbvector.h"
#include "cmbvectorset.h"
#include "cmbvectorsetset.h"
#include "MCMC.h"
 
enum class transformation { linear, softmax };
 
struct profiles_data
{
    vector<string> element_names;
    vector<string> sample_names;
    vector<vector<double>> values;
};
 
struct element_data
{
    string group_name;
    vector<string> sample_names;
    vector<double> values;
};
 
struct element_data_groups
{
    vector<string> group_names;
    vector<string> sample_names;
    vector<vector<double>> values;
};
 
struct DFA_result
{
    CMBVectorSet eigen_vectors;
    CMBVectorSet projected;
    CMBVector p_values;
    CMBVector wilkslambda;
    CMBVector F_test_P_value;
    CMBVectorSetSet multi_projected;
};
 
struct ANOVA_info
{
    double F = 0;
    double SST = 0;
    double SSB = 0;
    double SSW = 0;
    double MSB = 0;
    double MSW = 0;
    double p_value = 0;
};
 
enum class estimation_mode { only_contributions, elemental_profile_and_contribution, source_elemental_profiles_based_on_source_data };
 
enum class negative_reporting { list, complete };
 
class SourceSinkData : public map<string, Elemental_Profile_Set>
{
public:
    // ========== Construction and Assignment ==========
 
    SourceSinkData();
 
    SourceSinkData(const SourceSinkData& other);
 
    SourceSinkData& operator=(const SourceSinkData& other);
 
    // ========== Data Loading and Management ==========
 
    void Clear();
 
    Elemental_Profile_Set* AppendSampleSet(
        const string& name,
        const Elemental_Profile_Set& elemental_profile_set = Elemental_Profile_Set()
    );
 
    vector<string> SourceGroupNames() const;
 
    bool SetSelectedTargetSample(const string& sample_name);
 
    Elemental_Profile_Set* GetSampleSet(const string& name);
 
    bool SetTargetGroup(const string& targroup);
 
    string GetTargetGroup() const;
 
    // ========== Dataset Creation and Filtering ==========
 
    SourceSinkData CreateCorrectedDataset(
        const string& target,
        bool omnsizecorrect,
        map<string, element_information>* elementinfo
    );
 
    SourceSinkData CreateCorrectedAndFilteredDataset(
        bool exclude_samples,
        bool exclude_elements,
        bool omnsizecorrect,
        const string& target = ""
    ) const;
 
    SourceSinkData ExtractChemicalElements(bool isotopes) const;
 
    SourceSinkData BoxCoxTransformed(bool calculate_optimal_lambda = false);
 
    // ========== Element Information Management ==========
 
    map<string, element_information>* GetElementInformation();
 
    void IncludeExcludeElementsBasedOn(const vector<string>& elements);
 
    Distribution* GetFittedDistribution(const string& element_name);
    
    // ========== OM/Size Corrections ==========
 
    bool PerformRegressionVsOMAndSize(
        const string& om,
        const string& particle_size,
        regression_form form,
        const double& p_value_threshold = 0.05
    );
 
    vector<string> OMandSizeConstituents();
 
    vector<ResultItem> GetMLRResults();
 
    string FirstOMConstituent();
 
    string FirstSizeConstituent();
 
    void SetOMandSizeConstituents(const string& _omconstituent, const string& _sizeconsituent);
 
    void SetOMandSizeConstituents(const vector<string>& _omsizeconstituents);
 
    // ========== Outlier Analysis ==========
 
    void OutlierAnalysisForAll(const double& lower_threshold = -3, const double& upper_threshold = 3);
 
    // ========== Validation Tests ==========
 
    vector<string> NegativeValueCheck();
 
    CMBVector BracketTest(const string& target_sample, bool correct_based_on_om_n_size);
 
    CMBMatrix BracketTest(bool correct_based_on_om_n_size, bool exclude_elements, bool exclude_samples);
 
    // ========== Statistical Analysis - Discriminant Function Analysis ==========
 
    DFA_result DiscriminantFunctionAnalysis();
 
    DFA_result DiscriminantFunctionAnalysis(const string& source1, const string& source2);
 
    DFA_result DiscriminantFunctionAnalysis(const string& source1);
 
    vector<CMBVector> StepwiseDiscriminantFunctionAnalysis(
        const string& source1,
        const string& source2
    );
 
    vector<CMBVector> StepwiseDiscriminantFunctionAnalysis();
 
    vector<CMBVector> StepwiseDiscriminantFunctionAnalysis(const string& source1);
 
    // ========== Statistical Analysis - Differentiation Power ==========
 
    Elemental_Profile_Set DifferentiationPower(bool use_log, bool include_target);
 
    Elemental_Profile_Set DifferentiationPower_Percentage(bool include_target);
 
    Elemental_Profile_Set DifferentiationPower_P_value(bool include_target);
 
    Elemental_Profile DifferentiationPower_Percentage(const string& source1, const string& source2);
 
    Elemental_Profile DifferentiationPower(const string& source1, const string& source2, bool use_log);
 
    Elemental_Profile t_TestPValue(const string& source1, const string& source2, bool use_log);
 
    // ========== Statistical Analysis - ANOVA ==========
 
    CMBVector ANOVA(bool use_log);
 
    ANOVA_info ANOVA(const string& element, bool use_log);
 
    // ========== Parameter and Observation Setup ==========
 
    string GetParameterName(int index) const;
 
    bool InitializeParametersAndObservations(
        const string& targetsamplename,
        estimation_mode est_mode = estimation_mode::elemental_profile_and_contribution
    );
 
    void SetParameterEstimationMode(estimation_mode est_mode);
 
    bool SetParameterValue(size_t index, double value);
 
    void SetProgressReporter(ProgressReporter* _rtw);
 
    
 
    // ========== Optimization - Levenberg-Marquardt ==========
 
    bool SolveLevenberg_Marquardt(transformation trans = transformation::linear);
 
    CMBTimeSeriesSet LM_Batch(
        transformation transform,
        bool apply_om_size_correction,
        map<string, vector<string>>& negative_elements
    );
 
    CVector GetPredictedValues();
 
    CVector GradientUpdate(estimation_mode estmode = estimation_mode::elemental_profile_and_contribution);
 
    // ========== MCMC Analysis ==========
 
    Results MCMC(
        const string& target_sample,
        map<string, string> arguments,
        CMCMC<SourceSinkData>* mcmc,
        ProgressReporter* progress_window,
        const string& working_folder
    );
 
    CMBMatrix MCMC_Batch(
        map<string, string> arguments,
        CMCMC<SourceSinkData>* mcmc,
        ProgressReporter* progress_window,
        const string& working_folder,
        vector<string>* failed_writes = nullptr
    );
 
    // --- Likelihood and Objective Functions ---
 
    double GetObjectiveFunctionValue();
 
    double LogLikelihood(estimation_mode est_mode = estimation_mode::elemental_profile_and_contribution);
 
 
    // ========== Bootstrap and Validation ==========
 
    CMBTimeSeriesSet BootStrap(
        const double& percentage,
        unsigned int num_iterations,
        string target_sample,
        bool use_softmax
    );
 
    bool BootStrap(
        Results* results,
        const double& percentage,
        unsigned int num_iterations,
        string target_sample,
        bool use_softmax
    );
 
    CMBTimeSeriesSet VerifySource(
        const string& source_group,
        bool use_softmax,
        bool apply_om_size_correction
    );
 
    // ========== Results Retrieval ==========
 
    ResultItem GetContribution();
 
    ResultItem GetObservedElementalProfile();
 
    ResultItem GetObservedElementalProfile_Isotope();
 
    ResultItem GetPredictedElementalProfile(
        parameter_mode param_mode = parameter_mode::based_on_fitted_distribution
    );
 
    ResultItem GetPredictedElementalProfile_Isotope(
        parameter_mode param_mode = parameter_mode::based_on_fitted_distribution
    );
 
    ResultItem GetObservedvsModeledElementalProfile(
        parameter_mode param_mode = parameter_mode::based_on_fitted_distribution
    );
 
    ResultItem GetObservedvsModeledElementalProfile_Isotope(
        parameter_mode param_mode = parameter_mode::based_on_fitted_distribution
    );
 
    ResultItem GetCalculatedElementMeans();
 
    ResultItem GetCalculatedElementSigma();
 
    ResultItem GetCalculatedElementMu();
 
    ResultItem GetEstimatedElementMu();
 
    ResultItem GetEstimatedElementMean();
 
    ResultItem GetEstimatedElementSigma();
 
    vector<ResultItem> GetSourceProfiles();
 
    // ========== Utility ==========
 
    void AddtoToolsUsed(const string& tool);
 
    // ========== Accessors (inline for performance) ==========
 
    vector<Parameter>& Parameters();
 
    size_t ParametersCount();
 
    size_t ObservationsCount();
 
    Parameter* parameter(size_t i);
 
    const Parameter* parameter(size_t i) const;
 
    Observation* observation(size_t i);
 
    estimation_mode ParameterEstimationMode();
 
    ConcentrationSet* GetElementDistribution(const string& element_name);
 
    ConcentrationSet* GetElementDistribution(const string& element_name, const string& sample_group);
 
    element_information* GetElementInformation(const string& element_name);
 
    // --- Distribution Management ---
 
    void PopulateElementDistributions();
 
    void PopulateElementInformation(const map<string, element_information>* ElementInfo = nullptr);
 
    void AssignAllDistributions();
 
    // --- Data Extraction and Management ---
 
    int CountElements(bool exclude_elements) const;
 
    vector<string> GetSampleNames(const string& group_name) const;
 
    vector<string> GetGroupNames() const;
 
    vector<string> GetElementNames() const;
 
    void IncludeExcludeAllElements(bool include_in_analysis);
 
    Elemental_Profile Sample(const string& sample_name) const;
 
    profiles_data ExtractConcentrationData(const vector<vector<string>>& indicators) const;
 
    Elemental_Profile_Set ExtractSamplesAsProfileSet(const vector<vector<string>>& indicators) const;
 
    element_data ExtractElementConcentrations(const string& element, const string& group) const;
 
    map<string, vector<double>> ExtractElementDataByGroup(const string& element) const;
 
    SourceSinkData ExtractSpecificElements(const vector<string>& element_list) const;
 
    map<string, ConcentrationSet> ExtractConcentrationSet();
 
    // --- Element Ordering ---
 
    vector<string> ElementsToBeUsedInCMB();
 
    vector<string> IsotopesToBeUsedInCMB();
 
    void PopulateConstituentOrders();
 
    string SelectedTargetSample() const;
 
    // --- Tools Tracking ---
 
    bool ToolsUsed(const string& tool_name);
 
    // --- File I/O ---
 
    string GetOutputPath() const;
 
    QMap<QString, double>* GetOptions();
 
    bool SetOutputPath(const string& output_path);
 
    QJsonObject ElementInformationToJsonObject() const;
 
    QJsonArray ToolsUsedToJsonObject() const;
 
    QJsonObject OptionsToJsonObject() const;
 
    bool ReadToolsUsedFromJsonObject(const QJsonArray& jsonarray);
 
    bool ReadElementInformationfromJsonObject(const QJsonObject& jsonobject);
 
    bool ReadElementDatafromJsonObject(const QJsonObject& jsonobject);
 
    bool ReadOptionsfromJsonObject(const QJsonObject& jsonobject);
 
    QJsonObject ElementDataToJsonObject() const;
 
    bool WriteDataToFile(QFile* file);
 
    bool WriteToFile(QFile* file);
 
    bool ReadFromFile(QFile* fil);
private:
    // ========== Private Data Members ==========
 
    // Element Metadata and Distributions
    map<string, element_information> element_information_;
    map<string, ConcentrationSet> element_distributions_;
 
    // Counters
    int numberofconstituents_;
    int numberofisotopes_;
    int numberofsourcesamplesets_;
 
    // MCMC/Optimization Data
    vector<Observation> observations_;
    vector<Parameter> parameters_;
 
    // File Paths
    string outputpath_;
 
    // Group and Sample Identifiers
    string target_group_;
    string selected_target_sample_;
 
    // Ordering Vectors
    vector<string> samplesetsorder_;
    vector<string> constituent_order_;
    vector<string> element_order_;
    vector<string> isotope_order_;
    vector<string> size_om_order_;
 
    // Analysis Settings
    estimation_mode parameter_estimation_mode_;
    string omconstituent_;
    string sizeconsituent_;
    double regression_p_value_threshold_;
    double distance_coeff_;
    double error_stdev_;
    double error_stdev_isotope_;
    double epsilon_;
 
    // UI and Tracking
    ProgressReporter* rtw_ = nullptr;
    list<string> tools_used_;
    QMap<QString, double> options_;
 
    // ========== Private Helper Methods ==========
 
 
    Elemental_Profile* GetElementalProfile(const string& sample_name);
 
public:
    // --- Element Role Conversion ---
 
    QString Role(const element_information::role& role) const;
 
    element_information::role Role(const QString& role_string) const;
private:
 
    // --- Parameter Management ---
 
    double GetParameterValue(size_t index) const;
 
    bool SetParameterValue(const CVector& values);
 
    CVector GetParameterValue() const;
 
    // --- Contribution Initialization ---
 
    bool InitializeContributionsRandomly();
 
    bool InitializeContributionsRandomlySoftmax();
 
    // --- Contribution Accessors ---
 
    CVector GetSourceContributions();
 
    CVector GetContributionVector(bool include_all = true);
 
    CVector GetContributionVectorSoftmax();
 
    void SetContribution(size_t source_index, double contribution_value);
 
    void SetContributionSoftmax(size_t source_index, double softmax_value);
 
    void SetContribution(const CVector& contributions);
 
    void SetContributionSoftmax(const CVector& softmax_params);
 
    // --- Distribution Parameter Accessors ---
 
    Parameter* GetElementDistributionMuParameter(size_t element_index, size_t source_index);
 
    Parameter* GetElementDistributionSigmaParameter(size_t element_index, size_t source_index);
 
    double GetElementDistributionMuValue(size_t element_index, size_t source_index);
 
    double GetElementDistributionSigmaValue(size_t element_index, size_t source_index);
 
    // --- Observation Data ---
 
    CVector ObservedDataforSelectedSample(const string& SelectedTargetSample = "");
 
    CVector ObservedDataforSelectedSample_Isotope(const string& SelectedTargetSample = "");
 
    CVector ObservedDataforSelectedSample_Isotope_delta(const string& SelectedTargetSample = "");
 
    // --- Prediction Methods ---
 
    CVector PredictTarget(parameter_mode param_mode = parameter_mode::direct);
 
    CVector PredictTarget_Isotope(parameter_mode param_mode = parameter_mode::direct);
 
    CVector PredictTarget_Isotope_delta(parameter_mode param_mode = parameter_mode::based_on_fitted_distribution);
 
    CMatrix BuildSourceMeanMatrix(parameter_mode param_mode = parameter_mode::based_on_fitted_distribution);
 
    CMatrix BuildSourceMeanMatrix_Isotopes(parameter_mode param_mode = parameter_mode::based_on_fitted_distribution);
 
    // --- Likelihood and Prior Calculations ---
    double LogPriorContributions();
 
    double LogLikelihoodSourceElementalDistributions();
 
    double LogLikelihoodModelvsMeasured(estimation_mode est_mode = estimation_mode::elemental_profile_and_contribution);
 
    double LogLikelihoodModelvsMeasured_Isotope(estimation_mode est_mode = estimation_mode::elemental_profile_and_contribution);
 
    // --- Optimization Helper Methods ---
 
    CVector Gradient(const CVector& parameters, estimation_mode est_mode);
 
    CVector ResidualVector();
 
    CVector_arma ResidualVector_arma();
 
    CMatrix_arma ResidualJacobian_arma();
 
    CMatrix ResidualJacobian();
 
    CMatrix ResidualJacobian_softmax();
 
    CVector OneStepLevenberg_Marquardt(double lambda);
 
    CVector OneStepLevenberg_Marquardt_softmax(double lambda);
 
    vector<string> GetSourceOrder() const;
 
    vector<string> SamplesetsOrder();
 
    vector<string> ConstituentOrder();
 
    vector<string> ElementOrder();
 
    vector<string> IsotopeOrder();
 
    vector<string> SizeOMOrder();
 
    // --- Statistical Methods ---
 
    int TotalNumberofSourceSamples() const;
 
    double GrandMean(const string& element, bool use_log);
 
    Elemental_Profile_Set LumpAllProfileSets();
 
    // --- Box-Cox Transformation ---
 
    CMBVector OptimalBoxCoxParameters();
 
    // --- Bootstrap Helpers ---
 
    SourceSinkData RandomlyEliminateSourceSamples(const double& percentage);
 
    vector<string> AllSourceSampleNames() const;
 
    vector<string> RandomlypickSamples(const double& percentage) const;
 
    SourceSinkData ReplaceSourceAsTarget(const string& source_sample_name) const;
 
    // --- DFA Helper Methods ---
 
    CMatrix WithinGroupCovarianceMatrix();
 
    CMatrix BetweenGroupCovarianceMatrix();
 
    CMatrix TotalScatterMatrix();
 
    double WilksLambda();
 
    double DFA_P_Value();
 
    CMBVectorSet DFA_Projected();
 
    CMBVectorSet DFA_Projected(const string& source1, const string& source2);
 
    CMBVectorSet DFA_Projected(const string& source1, SourceSinkData* original);
 
    CMBVector DFA_eigvector();
 
    CMBVector DFA_weight_vector(const string& source1, const string& source2);
 
    CMBVector DFATransformed(const CMBVector& eigenvector, const string& source_group);
 
    Elemental_Profile_Set TheRest(const string& excluded_source);
 
    // --- Mean Calculation ---
 
    CMBVector MeanElementalContent(const string& group_name);
 
    CMBVector MeanElementalContent();
 
    
};
 
#endif // SOURCESINKDATA_H