GitHub PullRequestReviewCommentEvent at 2014-04-07 11:29:35
Repository: shogun
Actor: karlnapf
Payload: {"comment":{"url":"https://api.github.com/repos/shogun-toolbox/shogun/pulls/comments/11339649","id":11339649,"diff_hunk":"@@ -0,0 +1,394 @@\n+/*\n+ * Copyright (c) The Shogun Machine Learning Toolbox\n+ * Written (w) 2014 Wu Lin\n+ * All rights reserved.\n+ *\n+ * Redistribution and use in source and binary forms, with or without\n+ * modification, are permitted provided that the following conditions are met:\n+ *\n+ * 1. Redistributions of source code must retain the above copyright notice, this\n+ *    list of conditions and the following disclaimer.\n+ * 2. Redistributions in binary form must reproduce the above copyright notice,\n+ *    this list of conditions and the following disclaimer in the documentation\n+ *    and/or other materials provided with the distribution.\n+ *\n+ * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS \"AS IS\" AND\n+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED\n+ * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE\n+ * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR\n+ * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES\n+ * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;\n+ * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND\n+ * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT\n+ * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS\n+ * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n+ *\n+ * The views and conclusions contained in the software and documentation are those\n+ * of the authors and should not be interpreted as representing official policies,\n+ * either expressed or implied, of the Shogun Development Team.\n+ *\n+ */\n+\n+\n+#ifdef HAVE_EIGEN3\n+\n+#include <shogun/machine/gp/LogitPiecewiseBoundLikelihood.h>\n+\n+\n+using namespace shogun;\n+using namespace Eigen;\n+\n+CLogitPiecewiseBoundLikelihood::CLogitPiecewiseBoundLikelihood()\n+\t: CLogitLikelihood()\n+{\n+\tinit();\n+}\n+\n+CLogitPiecewiseBoundLikelihood::~CLogitPiecewiseBoundLikelihood()\n+{\n+}\n+\n+\n+void CLogitPiecewiseBoundLikelihood::set_bound(SGMatrix<float64_t> bound)\n+{\n+\tm_bound = bound;\n+}\n+\n+SGVector<float64_t> CLogitPiecewiseBoundLikelihood::get_variational_expection()\n+{\n+\tconst Map<VectorXd> eigen_c(m_bound.get_column_vector(0), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_b(m_bound.get_column_vector(1), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_a(m_bound.get_column_vector(2), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_mu(m_mu.vector, m_mu.vlen);\n+\tconst Map<VectorXd> eigen_s2(m_s2.vector, m_s2.vlen);\n+\tMap<VectorXd> eigen_l(m_bound.get_column_vector(3), m_bound.num_rows);\n+\tMap<VectorXd> eigen_h(m_bound.get_column_vector(4), m_bound.num_rows);\n+\n+\tindex_t num_rows = m_bound.num_rows; \n+\tindex_t num_cols = m_mu.vlen;\n+\n+\tconst Map<MatrixXd> eigen_cdf_diff(m_cdf_diff.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_pl(m_pl.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_ph(m_ph.matrix, num_rows, num_cols);\n+\n+\tfloat64_t l_bak = eigen_l(0);\n+\t//l(1) = 0; \n+\teigen_l(0) = 0;\n+\n+\tfloat64_t h_bak = eigen_h(eigen_h.size()-1);\n+\t//h(end) = 0; \n+\teigen_h(eigen_h.size()-1) = 0;\n+\n+\t//ex0 = ch-cl;\n+\tconst Map<MatrixXd> & eigen_ex0 = eigen_cdf_diff;\n+\n+\t//%ex1= v.*(pl-ph) + m.*(ch-cl);\n+\tMatrixXd eigen_ex1 = ((eigen_pl - eigen_ph).array().rowwise()*eigen_s2.array().transpose() \n+\t\t+ eigen_cdf_diff.array().rowwise()*eigen_mu.array().transpose()).matrix();\n+\n+\t//ex2 = bsxfun(@times, v, (bsxfun(@plus, l, m)).*pl - (bsxfun(@plus, h, m)).*ph) + bsxfun(@times, (v+m.^2), ex0);\n+\t\n+\t//bsxfun(@plus, l, m)).*pl - (bsxfun(@plus, h, m)).*ph\n+\tMatrixXd eigen_ex2 = ((eigen_mu.replicate(1,eigen_l.rows()).array().transpose().colwise() + eigen_l.array())*eigen_pl.array()\n+\t\t- (eigen_mu.replicate(1,eigen_h.rows()).array().transpose().colwise() + eigen_h.array())*eigen_ph.array()).matrix();\n+\n+\t//bsxfun(@times, v, (bsxfun(@plus, l, m).*pl - bsxfun(@plus, h, m).*ph\n+\teigen_ex2 = (eigen_ex2.array().rowwise()*eigen_s2.array().transpose()).matrix();\n+\n+\t//ex2 = bsxfun(@times, v, (bsxfun(@plus, l, m)).*pl - (bsxfun(@plus, h, m)).*ph) + bsxfun(@times, (v+m.^2), ex0);\n+\teigen_ex2 += (eigen_cdf_diff.array().rowwise()*(eigen_mu.array().pow(2)\n+\t\t\t+ eigen_s2.array()).transpose()).matrix();\n+\n+\tSGVector<float64_t> f(m_mu.vlen);\n+\tMap<VectorXd> eigen_f(f.vector, f.vlen);\n+\n+\t//%f = sum((a.*ex2 + b.*ex1 + c.*ex0),1);\n+\teigen_f = (eigen_ex2.array().colwise()*eigen_a.array() \n+\t\t+ eigen_ex1.array().colwise()*eigen_b.array() \n+\t\t+ eigen_ex0.array().colwise()*eigen_c.array()).colwise().sum().matrix();\n+\n+\n+\teigen_l(0) = l_bak;\n+\teigen_h(eigen_h.size()-1) = h_bak;\n+\n+\tMap<VectorXd>eigen_lab(m_lab.vector, m_lab.vlen);\n+\teigen_f = eigen_lab.cwiseProduct(eigen_mu) - eigen_f; \n+\treturn f;\n+}\n+\n+\n+SGVector<float64_t> CLogitPiecewiseBoundLikelihood::get_variational_first_derivative(\n+\t\tconst TParameter* param) const\n+{\n+\tREQUIRE(param, \"Param name is required (param should not be NULL)\\n\");\n+\tREQUIRE(!(strcmp(param->m_name, \"mu\") && strcmp(param->m_name, \"sigma2\")),\n+\t\t\"Can't compute derivative of the variational expection\", \n+\t\t\"of log LogitLikelihood using the piecewise bound\", \n+\t\t\"wrt %s.%s parameter\",\n+\t\tget_name(), param->m_name);\n+\n+\tconst Map<VectorXd> eigen_c(m_bound.get_column_vector(0), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_b(m_bound.get_column_vector(1), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_a(m_bound.get_column_vector(2), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_mu(m_mu.vector, m_mu.vlen);\n+\tconst Map<VectorXd> eigen_s2(m_s2.vector, m_s2.vlen);\n+\tMap<VectorXd> eigen_l(m_bound.get_column_vector(3), m_bound.num_rows);\n+\tMap<VectorXd> eigen_h(m_bound.get_column_vector(4), m_bound.num_rows);\n+\tindex_t num_rows = m_bound.num_rows; \n+\tindex_t num_cols = m_mu.vlen;\n+\tconst Map<MatrixXd> eigen_cdf_diff(m_cdf_diff.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_pl(m_pl.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_ph(m_ph.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_weighted_pdf_diff(m_weighted_pdf_diff.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_h2_plus_s2(m_h2_plus_s2.matrix, num_rows, num_cols);\n+\tconst Map<MatrixXd> eigen_l2_plus_s2(m_l2_plus_s2.matrix, num_rows, num_cols);\n+\tfloat64_t l_bak = eigen_l(0);\n+\t//l(1) = 0; \n+\teigen_l(0) = 0;\n+\tfloat64_t h_bak = eigen_h(eigen_h.size()-1);\n+\t//h(end) = 0; \n+\teigen_h(eigen_h.size()-1) = 0;\n+\n+\tSGVector<float64_t> result(m_mu.vlen);\n+\n+\tif (strcmp(param->m_name, \"mu\") ==0)\n+\t{\n+\t\t//bsxfun(@plus, bsxfun(@plus, l.^2, v), v).*pl - bsxfun(@plus, bsxfun(@plus, h.^2, v), v).*ph;\n+\t\tMatrixXd eigen_dmu2 = ((eigen_l2_plus_s2.array().rowwise()+eigen_s2.array().transpose())*eigen_pl.array()\n+\t\t\t- (eigen_h2_plus_s2.array().rowwise()+eigen_s2.array().transpose())*eigen_ph.array()).matrix();\n+\t\t//bsxfun(@times, ch-cl, (2.0*mu))\n+\t\teigen_dmu2 += (eigen_cdf_diff.array().rowwise()*(2.0*eigen_mu).array().transpose()).matrix();\n+\n+\t\tSGVector<float64_t> & gmu = result;\n+\t\tMap<VectorXd> eigen_gmu(gmu.vector, gmu.vlen);\n+\n+\t\t//gmu = bsxfun(@times, dmu2, _a) + bsxfun(@times, pl.*l - ph.*h + ch - chl, b) + bsxfun(@times, pl - ph, c)\n+\t\t//gmu = sum(gmu,1)\n+\t\teigen_gmu = ((eigen_dmu2.array().colwise()*eigen_a.array()) \n+\t\t\t+ ((eigen_weighted_pdf_diff + eigen_cdf_diff).array().colwise()*eigen_b.array())\n+\t\t\t+ ( (eigen_pl - eigen_ph).array().colwise()*eigen_c.array())).colwise().sum().matrix();\n+\n+\t\tMap<VectorXd>eigen_lab(m_lab.vector, m_lab.vlen);\n+\t\teigen_gmu = eigen_lab - eigen_gmu;\n+\t}\n+\telse\n+\t{\n+\n+\t\t//gv_0 = bsxfun(@plus, l, -mu).*pl - bsxfun(@plus, h, -mu).*ph;\n+\t\tMatrixXd eigen_gs2_0 = (((-eigen_mu).replicate(1,eigen_l.rows()).array().transpose().colwise() + eigen_l.array())*eigen_pl.array()\n+\t\t\t- ((-eigen_mu).replicate(1,eigen_h.rows()).array().transpose().colwise() + eigen_h.array())*eigen_ph.array()).matrix();\n+\t\t//gv_0 = bsxfun(times, gv_0, c);\n+\t\teigen_gs2_0 = (eigen_gs2_0.array().colwise()*eigen_c.array()).matrix();\n+\n+\t\t//tmpl = l2_plus_v - bsxfun(@times, l, mu)\n+\t\tMatrixXd tmpl = (eigen_l2_plus_s2 - (eigen_mu.replicate(1,eigen_l.rows()).array().transpose().colwise()*eigen_l.array()).matrix()\n+\t\t\t).cwiseProduct(eigen_pl);\n+\n+\t\t//tmph = h2_plus_v - bsxfun(@times, h, mu)\n+\t\tMatrixXd tmph = (eigen_h2_plus_s2 - (eigen_mu.replicate(1,eigen_h.rows()).array().transpose().colwise()*eigen_h.array()).matrix()\n+\t\t\t).cwiseProduct(eigen_ph);\n+\n+\t\t//gv_1 = bsxfun(@times, tmpl - tmph, b);\n+\t\tMatrixXd eigen_gs2_1 = ((tmpl - tmph).array().colwise()*eigen_b.array()).matrix();\n+\n+\t\t//gv_2 = bsxfun(@times, tmpl, l) - bsxfun(@times, tmph, h);\n+\t\tMatrixXd eigen_gs2_2 = (tmpl.array().colwise()*eigen_l.array() - tmph.array().colwise()*eigen_h.array()).matrix();\n+\n+\t\t//gv_2 = bsxfun(@times, gv_2, a);\n+\t\teigen_gs2_2 = (eigen_gs2_2.array().colwise()*eigen_a.array()).matrix();\n+\n+\t\tSGVector<float64_t> & gs2 = result;\n+\t\tMap<VectorXd> eigen_gs2(gs2.vector, gs2.vlen);\n+\n+\t\t//gv = (bsxfun(@times, ch - cl + 0.5*pl.*l - ph.*h, a) + bsxfun(@times, gv_0 + gv_1 + gv_2, 1.0/(2.0*v))\n+\t\t//gv = sum(gv,1);\n+\t\teigen_gs2 = ((eigen_cdf_diff + 0.5*eigen_weighted_pdf_diff).array().colwise()*eigen_a.array()\n+\t\t\t+ (eigen_gs2_0 + eigen_gs2_1 + eigen_gs2_2).array().rowwise()/(2.0*eigen_s2).array().transpose()\n+\t\t\t).colwise().sum().matrix();\n+\n+\t\t//gv = -gv\n+\t\teigen_gs2 = -eigen_gs2;\n+\t}\n+\teigen_l(0) = l_bak;\n+\teigen_h(eigen_h.size()-1) = h_bak;\n+\n+\n+\treturn result;\n+}\n+\n+void CLogitPiecewiseBoundLikelihood::set_distribution(SGVector<float64_t> mu,\n+\tSGVector<float64_t> s2, const CLabels* lab)\n+{\n+\tREQUIRE(lab, \"Labels are required (lab should not be NULL)\\n\");\n+\n+\tREQUIRE((mu.vlen==s2.vlen) && (mu.vlen==lab->get_num_labels()),\n+\t\t\"Length of the vector of means (%d), length of the vector of \"\n+\t\t\"variances (%d) and number of labels (%d) should be the same\\n\",\n+\t\tmu.vlen, s2.vlen, lab->get_num_labels());\n+\n+\tREQUIRE(lab->get_label_type()==LT_BINARY,\n+\t\t\"Labels must be type of CBinaryLabels\\n\");\n+\n+\tfor(index_t i = 0; i < s2.vlen; ++i)\n+\t\tASSERT(s2[i] > 0.0);\n+\n+\tm_lab = ((CBinaryLabels*)lab)->get_labels();\n+\t/** convert the input label to standard label used in the class\n+\t *\n+\t *  Note that Shogun uses  -1 and 1 as labels and this class uses \n+\t *  0 and 1 repectively.\n+\t *\n+\t */\n+\tfor(index_t i = 0; i < m_lab.size(); ++i)\n+\t\tm_lab[i] = CMath::max(m_lab[i], 0.0);\n+\n+\tm_mu = mu;\n+\n+\tm_s2 = s2;\n+\n+\tprecompute();\n+}\n+\n+void CLogitPiecewiseBoundLikelihood::init()\n+{\n+\tSG_ADD(&m_bound, \"bound\", \n+\t\t\"Variational piecewise bound for logit likelihood\",\n+\t\tMS_NOT_AVAILABLE);\n+\n+\tSG_ADD(&m_mu, \"mu\", \n+\t\t\"The mean of variational normal distribution\",\n+\t\tMS_AVAILABLE, GRADIENT_AVAILABLE);\n+\n+\tSG_ADD(&m_s2, \"sigma2\", \n+\t\t\"The variance of variational normal distribution\",\n+\t\tMS_AVAILABLE,GRADIENT_AVAILABLE);\n+\n+\tSG_ADD(&m_lab, \"y\", \n+\t\t\"The data/labels (must be 0 or 1) drawn from the distribution\",\n+\t\tMS_NOT_AVAILABLE);\n+\n+\tSG_ADD(&m_pl, \"pdf_l\", \n+\t\t\"The pdf given the lower range and parameters(mu and variance)\",\n+\t\tMS_NOT_AVAILABLE);\n+\tSG_ADD(&m_ph, \"pdf_h\", \n+\t\t\"The pdf given the higher range and parameters(mu and variance)\",\n+\t\tMS_NOT_AVAILABLE);\n+\tSG_ADD(&m_cdf_diff, \"cdf_h_minus_cdf_l\", \n+\t\t\"The CDF difference between the lower and higher range given the parameters(mu and variance)\",\n+\t\tMS_NOT_AVAILABLE);\n+\tSG_ADD(&m_l2_plus_s2, \"l2_plus_sigma2\", \n+\t\t\"The result of l^2 + sigma^2\",\n+\t\tMS_NOT_AVAILABLE);\n+\tSG_ADD(&m_h2_plus_s2, \"h2_plus_sigma2\", \n+\t\t\"The result of h^2 + sigma^2\",\n+\t\tMS_NOT_AVAILABLE);\n+\tSG_ADD(&m_weighted_pdf_diff, \"weighted_pdf_diff\", \n+\t\t\"The result of l*pdf(l_norm)-h*pdf(h_norm)\",\n+\t\tMS_NOT_AVAILABLE);\n+}\n+\n+void CLogitPiecewiseBoundLikelihood::precompute()\n+{\n+\tconst Map<VectorXd> eigen_c(m_bound.get_column_vector(0), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_b(m_bound.get_column_vector(1), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_a(m_bound.get_column_vector(2), m_bound.num_rows);\n+\tconst Map<VectorXd> eigen_mu(m_mu.vector, m_mu.vlen);\n+\tconst Map<VectorXd> eigen_s2(m_s2.vector, m_s2.vlen);\n+\tMap<VectorXd> eigen_l(m_bound.get_column_vector(3), m_bound.num_rows);\n+\tMap<VectorXd> eigen_h(m_bound.get_column_vector(4), m_bound.num_rows);\n+\n+\tindex_t num_rows = m_bound.num_rows; \n+\tindex_t num_cols = m_mu.vlen;\n+\n+\tm_pl = SGMatrix<float64_t>(num_rows,num_cols);\n+\tm_ph = SGMatrix<float64_t>(num_rows,num_cols);\n+\tm_cdf_diff = SGMatrix<float64_t>(num_rows,num_cols);\n+\tm_l2_plus_s2 = SGMatrix<float64_t>(num_rows,num_cols);\n+\tm_h2_plus_s2 = SGMatrix<float64_t>(num_rows,num_cols);\n+\tm_weighted_pdf_diff = SGMatrix<float64_t>(num_rows,num_cols);\n+\n+\tMap<MatrixXd> eigen_pl(m_pl.matrix, num_rows, num_cols);\n+\n+\tMap<MatrixXd> eigen_ph(m_ph.matrix, num_rows, num_cols);\n+\n+\tMap<MatrixXd> eigen_cdf_diff(m_cdf_diff.matrix, num_rows, num_cols);\n+\n+\tMap<MatrixXd> eigen_l2_plus_s2(m_l2_plus_s2.matrix, num_rows, num_cols);\n+\n+\tMap<MatrixXd> eigen_h2_plus_s2(m_h2_plus_s2.matrix, num_rows, num_cols);\n+\n+\tMap<MatrixXd> eigen_weighted_pdf_diff(m_weighted_pdf_diff.matrix, num_rows, num_cols);\n+\n+\n+\tSGMatrix<float64_t> zl(num_rows, num_cols);\n+\tMap<MatrixXd> eigen_zl(zl.matrix, num_rows, num_cols);\n+\n+\tSGMatrix<float64_t> zh(num_rows, num_cols);\n+\tMap<MatrixXd> eigen_zh(zh.matrix, num_rows, num_cols);\n+\n+\t//bsxfun(@minus,l,m)\n+\teigen_zl = ((-eigen_mu).replicate(1,eigen_l.rows()).array().transpose().colwise() + eigen_l.array()).matrix();\n+\t//bsxfun(@minus,h,m)\n+\teigen_zh = ((-eigen_mu).replicate(1,eigen_h.rows()).array().transpose().colwise() + eigen_h.array()).matrix();\n+\n+\tVectorXd eigen_s_inv = eigen_s2.array().sqrt().inverse().matrix(); ","path":"src/shogun/machine/gp/LogitPiecewiseBoundLikelihood.cpp","position":336,"original_position":334,"commit_id":"7a6114266ea656628055c3722548b06ffa046d21","original_commit_id":"7260af9ff511c12776d57ad24ea589f2f439cf4c","user":{"login":"karlnapf","id":699071,"avatar_url":"https://avatars.githubusercontent.com/u/699071?","gravatar_id":"5dd0cd1feee1624c281cdd1801b9ea81","url":"https://api.github.com/users/karlnapf","html_url":"https://github.com/karlnapf","followers_url":"https://api.github.com/users/karlnapf/followers","following_url":"https://api.github.com/users/karlnapf/following{/other_user}","gists_url":"https://api.github.com/users/karlnapf/gists{/gist_id}","starred_url":"https://api.github.com/users/karlnapf/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/karlnapf/subscriptions","organizations_url":"https://api.github.com/users/karlnapf/orgs","repos_url":"https://api.github.com/users/karlnapf/repos","events_url":"https://api.github.com/users/karlnapf/events{/privacy}","received_events_url":"https://api.github.com/users/karlnapf/received_events","type":"User","site_admin":false},"body":"ah I am sorry, of course","created_at":"2014-04-07T11:29:35Z","updated_at":"2014-04-07T11:29:35Z","html_url":"https://github.com/shogun-toolbox/shogun/pull/2106#discussion_r11339649","pull_request_url":"https://api.github.com/repos/shogun-toolbox/shogun/pulls/2106","_links":{"self":{"href":"https://api.github.com/repos/shogun-toolbox/shogun/pulls/comments/11339649"},"html":{"href":"https://github.com/shogun-toolbox/shogun/pull/2106#discussion_r11339649"},"pull_request":{"href":"https://api.github.com/repos/shogun-toolbox/shogun/pulls/2106"}}}}