reels.h

Go to the documentation of this file.

/* Mercury-Reels
 
    Copyright 2023 Banco Bilbao Vizcaya Argentaria, S.A.
 
    This product includes software developed at
 
    BBVA (https://www.bbva.com/)
 
      Licensed under the Apache License, Version 2.0 (the "License");
    you may not use this file except in compliance with the License.
    You may obtain a copy of the License at
 
    http://www.apache.org/licenses/LICENSE-2.0
 
      Unless required by applicable law or agreed to in writing, software
    distributed under the License is distributed on an "AS IS" BASIS,
    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    See the License for the specific language governing permissions and
    limitations under the License.
*/
#include <algorithm>
#include <map>
#include <math.h>
#include <set>
#include <string>
#include <string.h>
#include <time.h>
#include <vector>
#include <stdarg.h>
#include <cstdint>
 
#ifdef DEBUG
#include <iostream>
#include <fstream>
#include <sstream>
#endif
 
#ifdef TEST
#ifndef INCLUDED_CATCH2
#define INCLUDED_CATCH2
 
#include "catch.hpp"
 
#endif
#endif
 
 
#ifndef INCLUDED_REELS_TYPES
#define INCLUDED_REELS_TYPES
 
 
namespace reels
{
 
#define IMAGE_BUFF_SIZE         6136                    
#define PRIORITY_SEEN_FACTOR    2000000000              
#define DEFAULT_NUM_EVENTS      1000                    
#define MAX_SEQ_LEN_IN_PREDICT  1000                    
#define PREDICT_MAX_TIME        (100*365.25*24*3600)    
#define WEIGHT_PRECISION        10000                   
 
typedef uint64_t                        ElementHash;    
typedef std::string                     String;         
typedef const char *                    pChar;          
 
typedef time_t                          TimePoint;      
typedef struct tm                       TimeStruct;     
typedef double                          ExtFloat;       
 
 
struct ImageBlock {
    int size;                                           
    int block_num;                                      
 
    uint8_t buffer[IMAGE_BUFF_SIZE];                    
};
 
 
typedef std::vector<ImageBlock>         BinaryImage;    
typedef BinaryImage                    *pBinaryImage;   
 
 
struct BinEventPt {
    ElementHash e;                                      
    ElementHash d;                                      
    double w;                                           
 
    bool operator==(const BinEventPt &o) const {
        return e == o.e && d == o.d && round(WEIGHT_PRECISION*w) == round(WEIGHT_PRECISION*o.w);
    }
 
    bool operator<(const BinEventPt &o) const {
        return e < o.e || (e == o.e && d < o.d) || (e == o.e && d == o.d && round(WEIGHT_PRECISION*w) < round(WEIGHT_PRECISION*o.w));
    }
};
 
 
struct BinTransaction: BinEventPt {
    ElementHash c;                                      
    TimePoint t;                                        
};
 
 
struct EventStat {
    uint64_t seen;              
    uint64_t code;              
    uint64_t priority;          
};
 
 
typedef std::map<BinEventPt, EventStat> EventMap;
 
 
typedef std::map<uint64_t, BinEventPt> PriorityMap;
 
 
typedef std::map<uint64_t, uint64_t> EventCodeMap;
 
 
struct CodeInTreeStatistics {
    uint64_t n_succ_seen;       
    uint64_t n_succ_target;     
    uint64_t n_incl_seen;       
    uint64_t n_incl_target;     
    uint64_t sum_dep;           
    int      n_dep;             
};
 
 
typedef std::map<uint64_t, CodeInTreeStatistics> CodeInTreeStatMap;
 
 
struct StringUsage {
    uint64_t seen;  
    String str;     
};
 
 
typedef std::map<ElementHash, StringUsage> StringUsageMap;
 
 
typedef std::vector<ElementHash> ClientIDs;
 
 
typedef std::set<ElementHash> ClientIDSet;
 
 
typedef std::map<TimePoint, uint64_t> Clip;
 
 
typedef std::map<ElementHash, Clip> ClipMap;
typedef ClipMap * pClipMap; 
 
 
typedef std::map<ElementHash, TimePoint> TargetMap;
typedef TargetMap * pTargetMap;                 
 
 
typedef std::vector<double> TimesToTarget;
 
 
typedef std::set<uint64_t> CodeSet;
typedef CodeSet * pCodeSet;                 
 
 
struct OptimizeEvalItem {
    double      t_hat;      
    TimePoint   t_obs;      
    int         seq_len;    
};
 
 
bool compare_optimize_eval(const OptimizeEvalItem a, const OptimizeEvalItem b); 
 
 
typedef std::vector<OptimizeEvalItem>   OptimizeEval;
 
 
struct CodeScoreItem {
    uint64_t    code;   
    double      score;  
};
 
 
typedef std::vector <CodeScoreItem> CodeScores;
 
 
typedef std::map<uint64_t, int> ChildIndex;
 
 
struct CodeTreeNode {
    uint64_t   n_seen;      
    uint64_t   n_target;    
    ExtFloat   sum_time_d;  
    ChildIndex child;       
};
typedef CodeTreeNode * pCodeTreeNode;           
 
 
typedef std::vector<CodeTreeNode>   CodeTree;
typedef CodeTree * pCodeTree;                   
 
 
enum Transform {tr_undefined, tr_linear, tr_log};
 
 
enum Aggregate {ag_undefined, ag_mean, ag_minimax, ag_longest};
 
 
// Forward declaration of utilities used in other functions.
uint64_t MurmurHash64A (const void *key, int len);
bool image_put(pBinaryImage p_bi, void *p_data, int size);
bool image_get(pBinaryImage p_bi, int &c_block, int &c_ofs, void *p_data, int size);
 
 
class Logger {
 
    public:
 
        void log_printf(const char *fmt, ...) {
            va_list args;
            va_start(args, fmt);
            log_printf(fmt, args);
            va_end(args);
        }
 
        void log_printf(const char *fmt, va_list args) {
            char buffer[256];
 
            vsnprintf(buffer, sizeof(buffer), fmt, args);
 
            log = log + buffer;
        }
 
        String log = "";
};
 
 
class Clips;        // Forward declaration
 
class Events {
 
    public:
 
        Events() {}
 
        bool store_strings  = true;                     
        int  max_num_events = DEFAULT_NUM_EVENTS;       
 
 
        void insert_row(pChar p_e,
                        pChar p_d,
                        double w);
 
 
        bool define_event(pChar    p_e,
                          pChar    p_d,
                          double   w,
                          uint64_t code);
 
 
        String optimize_events (Clips &clips, TargetMap &targets, int num_steps = 10, int codes_per_step = 5, double threshold = 0.0001,
                                pCodeSet p_force_include = nullptr, pCodeSet p_force_exclude = nullptr, Transform x_form = tr_linear,
                                Aggregate agg = ag_longest, double p = 0.5, int depth = 1000, bool as_states = true,
                                double exp_decay = 0.00693, double lower_bound_p = 0.95, bool log_lift = true);
 
 
        bool score_model(double &score, double &targ_prop, CodeInTreeStatMap &codes_stat, bool calc_tree_stats, Clips &clips,
                         TargetMap &targets, EventCodeMap code_dict, Transform x_form, Aggregate agg, double p, int depth, bool as_states);
 
 
        CodeScores get_top_codes(CodeInTreeStatMap &codes_stat, double targ_prop, double exp_decay, double lower_bound_p,
                                 bool log_lift);
 
 
        double linear_correlation(OptimizeEval &ev) {
            ExtFloat s_h = 0, s_o = 0, sho = 0, ssh = 0, sso = 0;
            int n = 0;
 
            for (OptimizeEval::iterator it = ev.begin(); it != ev.end(); ++it) {
                if (it->t_obs != 0) {
                    s_h += it->t_hat;
                    s_o += it->t_obs;
                    sho += it->t_hat*it->t_obs;
                    ssh += it->t_hat*it->t_hat;
                    sso += it->t_obs*it->t_obs;
 
                    n++;
                }
            }
            if (n == 0)
                return 0;
 
            double d2 = (n*ssh - s_h*s_h)*(n*sso - s_o*s_o);
 
            if (d2 <= 1e-20)
                return 0;
 
            return (n*sho - s_h*s_o)/sqrt(d2);
        }
 
 
        bool load(pBinaryImage &p_bi);
 
 
        bool load(pBinaryImage &p_bi, int &c_block, int &c_ofs);
 
 
        bool save(pBinaryImage &p_bi);
 
 
        inline void set_max_num_events(int max_events) {
            max_num_events = max_events;
        }
 
 
        inline void set_store_strings(bool store) {
            store_strings = store;
        }
 
 
        inline ElementHash add_str(pChar p_str) {
            int ll = strlen(p_str);
 
            if (!ll)
                return 0;
 
            ElementHash hash = MurmurHash64A(p_str, ll);
 
            if (!store_strings)
                return hash;
 
            StringUsageMap::iterator it = names_map.find(hash);
 
            if (it != names_map.end())
                it->second.seen++;
            else {
                StringUsage su = {1, p_str};
 
                names_map[hash] = su;
            }
            return hash;
        }
 
 
        inline void erase_str(ElementHash hash) {
 
            if (store_strings) {
                StringUsageMap::iterator it = names_map.find(hash);
 
                if (it != names_map.end()) {
                    if (--it->second.seen == 0)
                        names_map.erase(it);
                }
            }
        }
 
 
        inline String get_str(ElementHash hash) {
 
            if (store_strings) {
                if (!hash)
                    return "";
 
                StringUsageMap::iterator it = names_map.find(hash);
 
                if (it != names_map.end())
                    return it->second.str;
            }
            return "\x04";
        }
 
 
        inline uint64_t event_code(BinEventPt &ept) {
            EventMap::iterator it = event.find(ept);
 
            if (it == event.end())
                return 0;
 
            return it->second.code;
        }
 
 
        inline int num_events() {
            return event.size();
        }
 
 
        inline EventMap::iterator events_begin() {
            return event.begin();
        }
 
 
        inline EventMap::iterator events_end() {
            return event.end();
        }
 
 
        inline EventMap::iterator events_next_after_find(BinEventPt &ept) {
            EventMap::iterator it = event.find(ept);
 
            if (it != event.end())
                ++it;
 
            return it;
        }
 
#ifndef TEST
    private:
#endif
 
        uint64_t priority_low = 0;
        uint64_t next_code    = 0;
 
        StringUsageMap names_map = {};
        EventMap       event     = {};
        PriorityMap    priority  = {};
};
 
 
class Clients {
 
    public:
 
        Clients() {}
 
 
        inline ElementHash hash_client_id(pChar p_cli) {
            int ll = strlen(p_cli);
 
            return ll == 0 ? 0 : MurmurHash64A(p_cli, ll);
        }
 
 
        void add_client_id(pChar p_cli);
 
 
        bool load(pBinaryImage &p_bi);
 
 
        bool load(pBinaryImage &p_bi, int &c_block, int &c_ofs);
 
 
        bool save(pBinaryImage &p_bi);
 
        ClientIDs   id     = {};    
        ClientIDSet id_set = {};    
};
 
 
class TimeUtil {
 
    public:
 
        TimeUtil() {}
 
        char time_format[128] = "%Y-%m-%d %H:%M:%S";    
 
 
        inline TimePoint get_time(pChar p_t) {
 
            TimeStruct ts = {0};
 
            if (strptime(p_t, time_format, &ts) == nullptr)
                return -1;
 
            return timegm(&ts);
        }
 
 
        inline void set_time_format(pChar fmt) {
            strncpy(time_format, fmt, sizeof(time_format) - 1);
        }
};
 
 
class Clips : public TimeUtil {
 
    public:
 
        Clips() {}
 
 
        Clips(Clients clients, Events events) : clients(clients), events(events) {}
 
 
        Clips(const ClipMap &clip_map) {
            clients = Clients();
            events  = Events();
 
            clips = clip_map;
        }
 
 
        Clips(Clips &o_clips) {
            pBinaryImage p_bi = new BinaryImage;
 
            o_clips.save(p_bi);
 
            load(p_bi);
 
            delete p_bi;
        }
 
 
        bool scan_event(pChar p_e,
                        pChar p_d,
                        double w,
                        pChar p_c,
                        pChar p_t);
 
 
        inline void insert_event(ElementHash client_hash,
                                 uint64_t    code,
                                 TimePoint   time_pt) {
 
            ClipMap::iterator it = clips.find(client_hash);
 
            if (it == clips.end()) {
                Clip clip = {};
 
                clip[time_pt] = code;
 
                clips[client_hash] = clip;
 
            } else
                it->second[time_pt] = code;
        };
 
 
        bool load(pBinaryImage &p_bi);
 
 
        bool save(pBinaryImage &p_bi);
 
 
        inline pClipMap clip_map() {
            return &clips;
        }
 
 
        inline uint64_t num_events() {
 
            uint64_t ret = 0;
 
            for (ClipMap::iterator it = clips.begin(); it != clips.end(); ++it)
                ret += it->second.size();
 
            return ret;
        }
 
        inline void collapse_to_states() {
            for (ClipMap::iterator it_client = clips.begin(); it_client != clips.end(); ++it_client) {
                uint64_t last_code = 0xA30BdefacedCabal;
                for (Clip::const_iterator it = it_client->second.cbegin(); it != it_client->second.cend();) {
                    uint64_t code = it->second;
                    if (code == last_code)
                        it_client->second.erase(it++);
                    else
                        ++it;
                    last_code = code;
                }
            }
        }
 
#ifndef TEST
    private:
#endif
 
        Clients clients;
        Events  events;
        ClipMap clips   = {};
};
 
 
class Targets   : public TimeUtil {
 
    public:
 
        Targets(pClipMap p_clips, TargetMap target) : p_clips(p_clips), target(target) {
            CodeTreeNode root = {0, 0, 0, {}};
            tree.push_back(root);
        }
 
 
        bool insert_target(pChar p_c, pChar p_t);
 
 
        bool fit(Transform x_form, Aggregate agg, double p, int depth, bool as_states);
 
 
        TimesToTarget predict();
 
 
        TimesToTarget predict(Clients clients);
 
 
        TimesToTarget predict(pClipMap p_clips);
 
 
        void verbose_predict_clip(const ElementHash &client,
                                  Clip              &clip,
                                  TimePoint         &obs_time,
                                  bool              &target_yn,
                                  int               &longest_seq,
                                  uint64_t          &n_visits,
                                  uint64_t          &n_targets,
                                  double            &targ_mean_t);
 
 
        bool load(pBinaryImage &p_bi);
 
 
        bool save(pBinaryImage &p_bi);
 
 
        inline int update_node(int idx_parent, uint64_t code, bool target, ExtFloat time_d) {
 
            if (idx_parent == 0) {      // The root node contains the prediction of the zero-length clip.
                tree[0].n_seen++;
                if (target) {
                    tree[0].n_target++;
                    tree[0].sum_time_d += time_d;
                }
            }
 
            ChildIndex::iterator it = tree[idx_parent].child.find(code);
 
            if (it != tree[idx_parent].child.end()) {
                int idx = it->second;
 
                tree[idx].n_seen++;
                if (target) {
                    tree[idx].n_target++;
                    tree[idx].sum_time_d += time_d;
                }
 
                return idx;
            }
 
            CodeTreeNode node = {1, target, time_d, {}};
 
            tree.push_back(node);
 
            int idx = tree.size() - 1;
 
            tree[idx_parent].child[code] = idx;
 
            return idx;
        }
 
 
        inline double normal_pdf(double x) {
            // https://stackoverflow.com/questions/10847007/using-the-gaussian-probability-density-function-in-c
 
            static const double inv_sqrt_2pi = 0.3989422804014327;
 
            return exp(-0.5*x*x)*inv_sqrt_2pi;
        }
 
 
        inline double normal_cdf(double x) {
            // https://stackoverflow.com/questions/2328258/cumulative-normal-distribution-function-in-c-c
            // https://cplusplus.com/reference/cmath/erfc/
 
            static const double m_sqrt_dot_5 = 0.7071067811865476;
 
            return 0.5*erfc(-x*m_sqrt_dot_5);
        }
 
 
        inline double agresti_coull_upper_bound(uint64_t n_hits, uint64_t n_total) {
            // https://github.com/msn0/agresti-coull-interval/blob/master/src/agresti.js
 
            double n_tilde = n_total + binomial_z_sqr;
            double p_tilde = (n_hits + binomial_z_sqr_div_2)/n_tilde;
            double a       = binomial_z*sqrt(p_tilde*(1 - p_tilde)/n_tilde);
 
            return p_tilde + a;
        }
 
 
        inline double agresti_coull_lower_bound(uint64_t n_hits, uint64_t n_total) {
            // https://github.com/msn0/agresti-coull-interval/blob/master/src/agresti.js
 
            double n_tilde = n_total + binomial_z_sqr;
            double p_tilde = (n_hits + binomial_z_sqr_div_2)/n_tilde;
            double a       = binomial_z*sqrt(p_tilde*(1 - p_tilde)/n_tilde);
 
            return p_tilde - a;
        }
 
 
        inline double predict_time(CodeTreeNode &node) {
 
            if (node.n_target <= 0)
                return PREDICT_MAX_TIME;
 
            double lb     = std::max(1e-4, agresti_coull_lower_bound(node.n_target, node.n_seen));
            double mu_hat = transform == tr_linear ? ((double) node.sum_time_d)/node.n_target : exp(((double) node.sum_time_d)/node.n_target);
 
            return mu_hat/lb;
        }
 
 
        inline double predict_clip(Clip clip) {
 
            int idx = 0, n = 0;
 
            double t[MAX_SEQ_LEN_IN_PREDICT];
 
            for (Clip::reverse_iterator it = clip.rbegin(); it != clip.rend(); it++) {
                ChildIndex::iterator jt = tree[idx].child.find(it->second);
 
                if (jt == tree[idx].child.end())
                    break;
 
                idx = jt->second;
 
                t[n++] = predict_time(tree[idx]);
            }
 
            if (n == 0)
                return predict_time(tree[0]);
 
            if (aggregate == ag_longest)
                return t[n - 1];
 
            double ret = t[0];
 
            if (aggregate == ag_mean) {
                for (int i = 1; i < n; i++)
                    ret += t[i];
 
                return ret/n;
            }
 
            for (int i = 1; i < n; i++)
                ret = std::min(ret, t[i]);
 
            return ret;
        }
 
 
        bool recurse_tree_stats(int depth, int idx, int parent_idx, uint64_t code, CodeInTreeStatMap &codes_stat);
 
 
        inline int num_targets() {
            return target.size();
        }
 
 
        inline int tree_size() {
            return tree.size();
        }
 
 
        inline pClipMap clip_map() {
            return p_clips;
        }
 
 
        inline pCodeTree p_tree() {
            return &tree;
        }
 
 
        inline pTargetMap p_target() {
            return &target;
        }
 
#ifndef TEST
    private:
#endif
 
        pClipMap   p_clips;
        TargetMap  target;
        CodeTree   tree                 = {};
        Transform  transform            = tr_undefined;
        Aggregate  aggregate            = ag_undefined;
        double     binomial_z           = 0;
        double     binomial_z_sqr       = 0;
        double     binomial_z_sqr_div_2 = 0;
        int        tree_depth           = 0;
};
 
} // namespace reels
 
#endif // ifndef INCLUDED_REELS_TYPES