トップ > R > Rの言語仕様とか色々参考にできそうな関数とかのメモ

Rの言語仕様とか色々参考にできそうな関数とかのメモ

R C

なんというS式

R Language Definitionのp3のシンボルオブジェクト。

> as.list(quote(x + y))
[[1]]
`+`

[[2]]
x

[[3]]
y

> as.list(quote((x + y) * z))
[[1]]
`*`

[[2]]
(x + y)

[[3]]
z

モデル式

言語オブジェクト

R 言語を構成する三つのオブジェクトの型がある。cal ls, expressions, そして names である。R は型 "expression"のオブジェクトを持つので、他の文脈中で言葉による表現式の使用を避けることにする。特に構文的に正しい表現は statements という名前で参照する。これらのオブジェクトは "call", "expression" そして"name" というモードをそれぞれ持つ。これらは表現式から quote 機構を用いて直接につくり出すことができ、as.list 関数と as.call 関数を用いてリストから、そしてリストへ、変換することができる。構文解析木の成分は標準的な添字操作を用いて抜き出すことができる。

R Language Definition

lmとかに渡すオブジェクトの中身がどうなっているのか、という話。

> mode(formula(y~x))
[1] "call"
> as.call(formula(y~x))
y ~ x
> as.list(formula(y~x))
[[1]]
`~`

[[2]]
y

[[3]]
x

> terms(formula(y~x))
y ~ x
attr(,"variables")
list(y, x)
attr(,"factors")
  x
y 0
x 1
attr(,"term.labels")
[1] "x"
attr(,"order")
[1] 1
attr(,"intercept")
[1] 1
attr(,"response")
[1] 1
attr(,".Environment")
<environment: R_GlobalEnv>
> str(terms(formula(y~x)))
Classes 'terms', 'formula' length 3 y ~ x
  ..- attr(*, "variables")= language list(y, x)
  ..- attr(*, "factors")= int [1:2, 1] 0 1
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "y" "x"
  .. .. ..$ : chr "x"
  ..- attr(*, "term.labels")= chr "x"
  ..- attr(*, "order")= int 1
  ..- attr(*, "intercept")= int 1
  ..- attr(*, "response")= int 1
  ..- attr(*, ".Environment")=<environment: R_GlobalEnv>

例

ここのやつがかなり簡単にできる。次数の高すぎる多項式フィッティングよくないという例を20行も行かないでplotまでやる。こんだけ簡単にやれるとちょっとびびる。

polynomial.lm <- function(x,y,q){
  xnam <- paste("I(x^", 1:q, ")",sep="")
  fmla <- as.formula(paste("y ~ ", paste(xnam, collapse= "+")))
  predict(lm(fmla))
}

n <- 10
x <- seq(from=0,to=1,length.out=n)
y <- sin(2*pi*x) + rnorm(n)

q <- 1:9

par(mfrow=c(3,3))
sapply(q,function(q){
  plot(x,y,main=paste("q = ",q,sep=""))
  lines(x,(function(x){sin(2*pi*x)})(x))
  lines(x,polynomial.lm(x,y,q),col="red")
})

f:id:syou6162:20161103155813p:plain

関数オブジェクト

R では関数はオブジェクトであり、他のオブジェクトと良く似た方法で操作できる。関数は三つの基本成分、形式的引数のリスト、本体、そして環境を持つ。引数リストはコンマで区切られた引数のリストである。

R Language Definition

クロージャーオブジェクトの三つの部分を formals, body そして environment を用いて取り出し、操作することが可能である (また三つのどれもが代入式の左辺項に使うことができる)。最後の一つはしばしば不必要な環境取り込み要素を取り除くのに使われる。as.list と as.function を用い、関数をリスト構造へ、そしてリスト構造から関数へ、変換する機能がある。これらは S との互換性を維持するために提供されているが、使わないことを勧める。

R Language Definition

定義

ふむ、本当に3つの要素から成っているようだ。

struct closxp_struct {
    struct SEXPREC *formals;
    struct SEXPREC *body;
    struct SEXPREC *env;
};

環境

環境は二つのものからなると考えることができる。フレーム、つまりシンボルと値の対の集りと、それを取り巻く環境とである。あるシンボルの値が必要になると、フレームが検査され、マッチするシンボルが見付かると、それが返される。さもなければその外側の環境が次に調べられる、等が繰り返される。環境は木構造を作り、一つの環境は複数の子環境を持つことができるが、親は唯一つだけである。

R Language Definition

言語自体のプログラミング

R は、少なくとも数学公式や C 風の制御構造になれた者に取って、Lisp よりはよりとっつき易いプログラミングへのインタフェイスを提供するが、そのエンジンは実際はかなり Lisp 風である。R は構文解析された表現式や関数への直接のアクセスを許し、それらを変更した後実行したり、さらにはまったく新しい関数を一からつくり出すことを許す。

R Language Definition

表現式の導関数を計算したり、係数ベクトルから多項式を生成したりといった、この機能の幾つかの標準的応用がある。しかしながら、また R のインタープリタ部分の動作にとりより基本的な使用法が存在する。これらの幾つかは、幾つかのモデリング・プロットルーチン中でつくり出されるmodel.frame への呼び出し (決して手際良いとはいえないかも知れないが) に於けるように、関数を他の関数の一部分として再使用するために本質的である。

R Language Definition

treeApply.c

static SEXP apply_call(SEXP object, SEXP expr, SEXP symbol, SEXP env)
{
    if(TYPEOF(expr) == LANGSXP &&
       /* unfortunately, R uses the contents of an object to detect
	  missing arguments, meaning that assigning R_MissingArg to
	  anything totally confuses any function calls:  */
	object != R_MissingArg) {
	defineVar(symbol, object, env);
	return eval(expr, env);
    }
    else return object;
}

Rinternals.h

/* Closure Access Macros */
#define FORMALS(x)	((x)->u.closxp.formals)
#define BODY(x)		((x)->u.closxp.body)
#define CLOENV(x)	((x)->u.closxp.env)
#define DEBUG(x)	((x)->sxpinfo.debug)
#define SET_DEBUG(x,v)	(((x)->sxpinfo.debug)=(v))

typedef struct SEXPREC {
    SEXPREC_HEADER;
    union {
	struct primsxp_struct primsxp;
	struct symsxp_struct symsxp;
	struct listsxp_struct listsxp;
	struct envsxp_struct envsxp;
	struct closxp_struct closxp;
	struct promsxp_struct promsxp;
    } u;
} SEXPREC, *SEXP;

struct closxp_struct {
    struct SEXPREC *formals;
    struct SEXPREC *body;
    struct SEXPREC *env;
};

struct listsxp_struct {
    struct SEXPREC *carval;
    struct SEXPREC *cdrval;
    struct SEXPREC *tagval;
};

kmeans.c

Rでの関数定義はこんな感じ。

kmeans <-
function(x, centers, iter.max = 10, nstart = 1,
         algorithm = c("Hartigan-Wong", "Lloyd", "Forgy", "MacQueen"))

引数の意味はそれぞれこんな感じ。

Arguments:

       x: A numeric matrix of data, or an object that can be coerced to
          such a matrix (such as a numeric vector or a data frame with
          all numeric columns).

 centers: Either the number of clusters or a set of initial (distinct)
          cluster centres.  If a number, a random set of (distinct)
          rows in 'x' is chosen as the initial centres.

iter.max: The maximum number of iterations allowed.

  nstart: If 'centers' is a number, how many random sets should be
          chosen?

algorithm: character: may be abbreviated.

kmeansのexample。

     x <- rbind(matrix(rnorm(100, sd = 0.3), ncol = 2),
                matrix(rnorm(100, mean = 1, sd = 0.3), ncol = 2))
     colnames(x) <- c("x", "y")
     (cl <- kmeans(x, 2))

kmeans.Rの途中はこんな感じになっている。

                       Z <- .C("kmeans_Lloyd", as.double(x), as.integer(m),
                               as.integer(ncol(x)),
                               centers = as.double(centers), as.integer(k),
                               c1 = integer(m), iter = as.integer(iter.max),
                               nc = integer(k), wss = double(k),
                               PACKAGE="stats")
                       if(Z$iter > iter.max)
                           warning("did not converge in ",
                                  iter.max, " iterations", call.=FALSE)
                       if(any(Z$nc == 0))
                           warning("empty cluster: try a better set of initial centers", call.=FALSE)
                       Z

kmeans_Lloydの定義。当然パラメータは全てポインタで渡す。

void kmeans_Lloyd(double *x, int *pn, int *pp, double *cen, int *pk, int *cl,
		  int *pmaxiter, int *nc, double *wss)
{
    int n = *pn, k = *pk, p = *pp, maxiter = *pmaxiter;
    int iter, i, j, c, it, inew = 0;
    double best, dd, tmp;
    Rboolean updated;

    for(i = 0; i < n; i++) cl[i] = -1;
    for(iter = 0; iter < maxiter; iter++) {
	updated = FALSE;
	for(i = 0; i < n; i++) {
	    /* find nearest centre for each point */
	    best = R_PosInf;
	    for(j = 0; j < k; j++) {
		dd = 0.0;
		for(c = 0; c < p; c++) {
		    tmp = x[i+n*c] - cen[j+k*c];
		    dd += tmp * tmp;
		}
		if(dd < best) {
		    best = dd;
		    inew = j+1;
		}
	    }
	    if(cl[i] != inew) {
		updated = TRUE;
		cl[i] = inew;
	    }
	}
	if(!updated) break;
	/* update each centre */
	for(j = 0; j < k*p; j++) cen[j] = 0.0;
	for(j = 0; j < k; j++) nc[j] = 0;
	for(i = 0; i < n; i++) {
	    it = cl[i] - 1; nc[it]++;
	    for(c = 0; c < p; c++) cen[it+c*k] += x[i+c*n];
	}
	for(j = 0; j < k*p; j++) cen[j] /= nc[j % k];
    }

    *pmaxiter = iter + 1;
    for(j = 0; j < k; j++) wss[j] = 0.0;
    for(i = 0; i < n; i++) {
	it = cl[i] - 1;
	for(c = 0; c < p; c++) {
	    tmp = x[i+n*c] - cen[it+k*c];
	    wss[it] += tmp * tmp;
	}
    }
}

スカラーな(?)ポインタは最初に全部ポインタじゃない変数に変換してある。

    int n = *pn, k = *pk, p = *pp, maxiter = *pmaxiter;

Rでの行列をCではどう扱うか。2次元配列として、ではなく1次元の配列で頑張らないといけないことが見てとれる。

		    tmp = x[i+n*c] - cen[j+k*c];

Learning R

Learning R