Rcpp for everyone

3DQQGPSFWFSZPOF ௡ବ@teuder

"CPVUNF ௡ాਅथ .BTBLJ&5TVEB 1I% -JGFTDJFODFT #VTJOFTT$POTVMUBOU XJUI %BUB4DJFODF

3DQQ 3DQQJTBQBDLBHFUIBUFOBCMFTZPV UPJNQMFNFOU3GVODUJPOTJO$ #ZVTJOH3DQQ ZPVDBOXSJUF$ DPEFJOBTUZMFTJNJMBSUP3

%JSL&EEFMCVFUUFM

5PEBZ *EPOPUUBMLBCPVUIPXUPVTF3DQQ CVUJOUSPEVDFZPVTPNFVTFGVMSFTPVSDFT XIFOZPVXBOUUPMFBSO3DQQ

5IFNBKPSSFTPVSDFTGPS MFBSOJOH3DQQ "EWBODFE3 3DQQ(BMMFSZ 3DQQ/PUF 7JHOFUUFT

"EWBODFE3 (SFBUJOUSPEVDUJPOUP3DQQ #VUJUMBDLTTPNFEFUBJMT )PXDBO*NPWFUPUIFOFYUTUFQ

3DQQ(BMMFSZ (SFBUDPNQJMBUJPOPGHPPEFYBNQMFT #VUUIFJOGPSNBUJPOJTGSBHNFOUFE )PXDBO*HFUUJEZLOPXMFEHFBCPVU3DQQ

3DQQ/PUF (SFBU$ TUZMFSFGFSFODFEPDVNFOU #VUJUPOMZQSPWJEFMJNJUFEJOGPSNBUJPO )PXDBO*HFUNPSFQSBDUJDBM BOENPSFTZTUFNBUJDBMLOPXMFEHF

*OUIFTFEBZT *DPVMEOPUpOEBEPDVNFOUMJLFUIBU*XBOU #VUPOFEBZ *SFNFNCFSFEBQISBTF *`WFTFFOPOUIFJOUFSOFU

If there does not exist what you want, just create
it. That should have been the oath of HOXO-M ! ͳ͚Ε͹࡞Δɻ ͦΕ͕ϗΫιΤϜͷ੤͍ͩͬͨ͸ͣʂ

IUUQTXXXHJUCPPLDPNCPPLUFVEFSJOUSPEVDUJPOUPSDQQEFUBJMT *XSPUFB3DQQ8FC#PPL

&OHMJTIWFSTJPOJTBMTPBWBJMBCMF IUUQTXXXHJUCPPLDPNCPPLUFVEFSJOUSPEVDUJPOUPSDQQEFUBJMT

$POUFOUTPGUIFCPPL 7FDUPS .BUSJY -JTU %BUB'SBNF 4 3MJLFGVODUJPOT 1SPCBCJMJUZEJTUSJCVUJPOT 6TJOH3GVODUJPOTGSPN3DQQ "UUSJCVUFT
/"/BO*OG &OWJSPONFOU &SSPSIBOEMJOH BOENPSFʜ

4JODFUIF&OHMJTIWFSTJPO XBTSFMFBTFE 1BHFWJFXTEPVCMFEBOEJTTUJMMHSPXJOH 3FMFBTFPG&OHWFS /PW 5JNFTFSJFTPGQBHFWJFXTJOUIFMBTUNPOUIT

1BHFWJFXTCZDPVOUSJFT (FSNBOZJTUIFSEQMBDF #VUMBSHFMZCFIJOEGSPN64BOE+BQBO 1BHFWJFXTJOUIFMBTUNPOUIT

8IBUXPVMECF UIFHPPETJUVBUJPOUPVTF3DQQ w -PPQTJOXIJDIMBUFSJUFSBUJPOEFQFOETPOUIF QSFWJPVTJUFSBUJPO w 3FDVSSFOUDBMMTPGGVODUJPOJOMPPQ w "DDFTTJOHFBDIFMFNFOUTPGWFDUPSEBUBGSBNF w
$IBOHJOHTJ[FPGBWFDUPSBUSVOUJNF w 0QFSBUJPOUIBUOFFEBEWBODFEEBUBTUSVDUVSF BOEBMHPSJUIN

#FDBVTF3DQQSFRVJSFTZPVBOFYUSBDPTU :PVNJHIUCFBCMFUPNBYJNJ[FZPVSRVBMJUZPGMJGF XIFOZPVVTFJUUPCSFBLCPUUMFOFDLPGZPVSDPEF 8IBUXPVMECF UIFHPPETJUVBUJPOUPVTF3DQQ

*ONZDBTF $POWFSUJOHEBUBGSBNFUPTQBSTFNBUSJY df %>% as.matrix %>% Matrix::Matrix(sparse = TRUE) #VUDPOWFSUJOHNBUSJYGBJMTXIFOUIFEBUBJT
SFMBUJWFMZMBSHF BCPVUTFWFSBMUFOTPG(# df %>% asSparseMatrix 4P*DSFBUFB3DQQGVODUJPODPOWFSUJOHUIFEBUB EJSFDUMZUPTQBSTFNBUSJY

#include <Rcpp.h> using namespace Rcpp; %FpOJOHBGVODUJPO // [[Rcpp::export]] S4 asSparseMatrix(
DataFrame df ){ // prerequisite : // all the elements DataFrame is numeric/integer // and not containing NAs. // number of rows and columns int nrow = df.nrows(); int ncol = df.length();

6TJOHTUEWFDUPSJOTUFBEPG3DQQ7FDUPS #FDBVTFDIBOHJOHTJ[FPGWFDUPSBUSVOUJNF JTOPUFGpDJFOUJO3DQQ7FDUPS std::vector<R_xlen_t> rows; std::vector<R_xlen_t> cols; std::vector<double> vals;
4UPSJOHUIFQPTJUJPOTBOEWBMVFTPG OPO[FSPFMFNFOUTJOUIF%BUB'SBNF for(R_xlen_t col=0; col<ncol; ++col){ NumericVector column = df[col]; for(R_xlen_t row=0; row<nrow; ++row){ if(column[row]!=0.0){ rows.push_back(row+1); cols.push_back(col+1); vals.push_back(column[row]); } } }

$BMMJOH.BUSJYTQBSTF.BUSJY Environment env = Environment::namespace_env("Matrix"); Function sparseMatrix = env["sparseMatrix"];
$POWFSUJOHTUEWFDUPSUP/VNFSJD7FDUPS S4 sm = sparseMatrix( Named("i") = wrap(rows), Named("j") = wrap(cols), Named("x") = wrap(vals), Named("dims") = NumericVector::create(nrow,ncol)); 4FUUJOHSPXOBNFTBOEDPMOBNFT List dimnames = List::create(R_NilValue, df.names()); sm.attr("Dimnames") = dimnames; 3FUVSOJOHUIFTQBSTFNBUSJY return sm; }

$PODMVTJPO w :PVDBOJNQMFNFOU3GVODUJPOJO$ XJUI3DQQ w 5IFSFBSFTV⒏DJFOUSFTPVSDFTGPSMFBSOJOH3DQQ w 3DQQJTTVJUBCMFUPCSFBLCPUUMFOFDLTPGZPVSDPEF *IPQFNZl*OUSPEVDUJPOUP3DQQz XPVMEBJECSPBESBOHFPGQFPQMF
&OKPZ

"EWFSUJTJOH 605 1 2 3 4 5 6 7 8
9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 24 ষɹɹRcpp Rcpp Rcpp͸Rͷؔ਺ΛC++Ͱ࣮૷Ͱ͖ΔύοέʔδͰ͢ɻRͱྨࣅͨ͠ελΠϧͰهड़Ͱ͖ ΔΑ͏ʹ࣮૷͞Ε͍ͯΔͨΊɺC++ʹਂ͍஌͕ࣝͳͯ͘΋ར༻͠΍͘͢ͳ͍ͬͯ·͢ɻ͠ ͔΋ɺͦͷͨΊͷ࣮ߦ଎౓͸٘ਜ਼ʹ͞Ε͍ͯͳ͍ͷͰɺ୭Ͱ΋ϋΠύϑΥʔϚϯεͳ݁ՌΛ ಘΔ͜ͱ͕Ͱ͖·͢ɻ 24-1 Rcppͷ׆༻γʔϯ ࣍ͷΑ͏ͳέʔε͸C++Ͱ࣮૷͢Δ͜ͱʹΑΓɺRͱൺ΂ͯߴ଎Խ͕ݟࠐΊ·͢ɻ ɾ ܁Γฦ͠ॲཧɺಛʹ࣍ͷॲཧ͕લͷॲཧʹґଘ͓ͯ͠ΓฒྻԽͰ͖ͳ͍ ɾ ϕΫ τϧ΍ߦྻͷݸʑͷཁૉ΁ΞΫηε͢Δඞཁ͕͋Δ ɾ ϕΫ τϧͷαΠζΛಈతʹมߋ͍ͨ͠ ɾ ߴ౓ͳσʔλߏ଄΍ΞϧΰϦζϜΛ༻͍ͨॲཧΛߦ͍͍ͨ RcppͷύϑΥʔϚϯεΛࣔͨ͢Ίɺ܁Γฦ͠ॲཧͷྫͱͯ͠MCMCΞϧΰϦζϜͷҰछͰ͋ ΔΪϒεαϯϓϥʔ ʢ஫1ʣ ͷ࣮૷ྫΛࣔ͠·͢ɻ͜ͷྫͰ͸ɺΪϒεαϯϓϥʔʹΑΓඪ४2ม਺ ਖ਼ن෼෍͔Βn఺αϯϓϦϯά͍ͯ͠·͢ɻ ·ͣ͸ൺֱͷͨΊRΛ༻͍࣮ͨ૷ྫΛࣔ͠·͢ ʢϦετ24.1ʣ ɻ Ϧε τ24.1ɹGibbs.R GibbsR <- function(b, n, t){ # 2ม਺ඪ४ਖ਼ن෼෍͔Βn఺αϯϓϦ ϯά # b : 2ม਺ͷڞ෼ࢄ # n : αϯϓϧ਺ # t : αϯϓϦ ϯάࣺͤͣͯΔִؒ X <- matrix(0, nrow = n, ncol = 2) x1 <- x2 <- 0 sd <- sqrt(1-b^2) for(i in 1:n){ for(j in 1:t){ x1 <- rnorm(1, b*x2, sd) x2 <- rnorm(1, b*x1, sd) ʢ஫1ʣ ߴ࣍ݩͷෳࡶͳ֬཰෼෍ʹै͏ཚ਺Λੜ੒͢ΔϚϧίϑ࿈࠯ϞϯςΧϧϩ๏ͱݺ͹ΕΔΞϧΰϦζϜͷҰछͰ͢ɻ 24-1 24 ষ ٕज़ධ࿦ࣾɹ੫ࠐԁ

Rcpp for everyone

Rcpp for everyone

teuder

More Decks by teuder

Other Decks in Technology

Featured

Transcript

3DQQGPSFWFSZPOF ௡ବ@teuder

"CPVUNF ௡ాਅथ .BTBLJ&5TVEB 1I% -JGFTDJFODFT #VTJOFTT$POTVMUBOU XJUI %BUB4DJFODF

3DQQ 3DQQJTBQBDLBHFUIBUFOBCMFTZPV UPJNQMFNFOU3GVODUJPOTJO$ #ZVTJOH3DQQ ZPVDBOXSJUF$ DPEFJOBTUZMFTJNJMBSUP3

%JSL&EEFMCVFUUFM

5PEBZ *EPOPUUBMLBCPVUIPXUPVTF3DQQ CVUJOUSPEVDFZPVTPNFVTFGVMSFTPVSDFT XIFOZPVXBOUUPMFBSO3DQQ

5IFNBKPSSFTPVSDFTGPS MFBSOJOH3DQQ "EWBODFE3 3DQQ(BMMFSZ 3DQQ/PUF 7JHOFUUFT

"EWBODFE3 (SFBUJOUSPEVDUJPOUP3DQQ #VUJUMBDLTTPNFEFUBJMT )PXDBO*NPWFUPUIFOFYUTUFQ

3DQQ(BMMFSZ (SFBUDPNQJMBUJPOPGHPPEFYBNQMFT #VUUIFJOGPSNBUJPOJTGSBHNFOUFE )PXDBO*HFUUJEZLOPXMFEHFBCPVU3DQQ

3DQQ/PUF (SFBU$ TUZMFSFGFSFODFEPDVNFOU #VUJUPOMZQSPWJEFMJNJUFEJOGPSNBUJPO )PXDBO*HFUNPSFQSBDUJDBM BOENPSFTZTUFNBUJDBMLOPXMFEHF

OUIFTFEBZT DPVMEOPUpOEBEPDVNFOUMJLFUIBUXBOU #VUPOFEBZ SFNFNCFSFEBQISBTF *`WFTFFOPOUIFJOUFSOFU

If there does not exist what you want, just create

IUUQTXXXHJUCPPLDPNCPPLUFVEFSJOUSPEVDUJPOUPSDQQEFUBJMT *XSPUFB3DQQ8FC#PPL

&OHMJTIWFSTJPOJTBMTPBWBJMBCMF IUUQTXXXHJUCPPLDPNCPPLUFVEFSJOUSPEVDUJPOUPSDQQEFUBJMT

$POUFOUTPGUIFCPPL 7FDUPS .BUSJY -JTU %BUB'SBNF 4 3MJLFGVODUJPOT 1SPCBCJMJUZEJTUSJCVUJPOT 6TJOH3GVODUJPOTGSPN3DQQ "UUSJCVUFT

4JODFUIF&OHMJTIWFSTJPO XBTSFMFBTFE 1BHFWJFXTEPVCMFEBOEJTTUJMMHSPXJOH 3FMFBTFPG&OHWFS /PW 5JNFTFSJFTPGQBHFWJFXTJOUIFMBTUNPOUIT

1BHFWJFXTCZDPVOUSJFT (FSNBOZJTUIFSEQMBDF #VUMBSHFMZCFIJOEGSPN64BOE+BQBO 1BHFWJFXTJOUIFMBTUNPOUIT

8IBUXPVMECF UIFHPPETJUVBUJPOUPVTF3DQQ w -PPQTJOXIJDIMBUFSJUFSBUJPOEFQFOETPOUIF QSFWJPVTJUFSBUJPO w 3FDVSSFOUDBMMTPGGVODUJPOJOMPPQ w "DDFTTJOHFBDIFMFNFOUTPGWFDUPSEBUBGSBNF w

#FDBVTF3DQQSFRVJSFTZPVBOFYUSBDPTU :PVNJHIUCFBCMFUPNBYJNJ[FZPVSRVBMJUZPGMJGF XIFOZPVVTFJUUPCSFBLCPUUMFOFDLPGZPVSDPEF 8IBUXPVMECF UIFHPPETJUVBUJPOUPVTF3DQQ

*ONZDBTF $POWFSUJOHEBUBGSBNFUPTQBSTFNBUSJY df %>% as.matrix %>% Matrix::Matrix(sparse = TRUE) #VUDPOWFSUJOHNBUSJYGBJMTXIFOUIFEBUBJT

#include <Rcpp.h> using namespace Rcpp; %FpOJOHBGVODUJPO // [[Rcpp::export]] S4 asSparseMatrix(

6TJOHTUEWFDUPSJOTUFBEPG3DQQ7FDUPS #FDBVTFDIBOHJOHTJ[FPGWFDUPSBUSVOUJNF JTOPUFGpDJFOUJO3DQQ7FDUPS std::vector<R_xlen_t> rows; std::vector<R_xlen_t> cols; std::vector<double> vals;

$BMMJOH.BUSJYTQBSTF.BUSJY Environment env = Environment::namespace_env("Matrix"); Function sparseMatrix = env["sparseMatrix"];

$PODMVTJPO w :PVDBOJNQMFNFOU3GVODUJPOJO$ XJUI3DQQ w 5IFSFBSFTV⒏DJFOUSFTPVSDFTGPSMFBSOJOH3DQQ w 3DQQJTTVJUBCMFUPCSFBLCPUUMFOFDLTPGZPVSDPEF IPQFNZlOUSPEVDUJPOUP3DQQz XPVMEBJECSPBESBOHFPGQFPQMF

"EWFSUJTJOH 605 1 2 3 4 5 6 7 8