Spaces and Convergences

Metric spaces: a set where a notion of distance (called a metric) between elements of the set is defined.

Vector spaces or Linear space: a mathematical structure formed by a collection of vectors.

Normed vector spaces: a Vector space with a norm defined. The norm is an abstraction of our usual concept of length.

Any normed vector space is a metric space by defining d(x,y)=||y-x||.

Metric space may not be a vector space.

Inner product spaces: Vector spaces endowed with inner product.

An inner product naturally induces an associated norm, thus an inner product space is also a normed vector space.

A normed vector spaces is also an inner product spaces if and only if its norm satisfies the parallelogram law.

Banach spaces: complete normed vector spaces.

Hilbert spaces: complete inner product spaces.

Weak convergence(in normed vector space)

Suppose X is a normed vector space, X^{*} is the continuous dual of X , and x_{1},x_{2},\cdots is a sequence in X, x_{o}\in X. Then we say that {x_{n}} converges weakly to x_{0} if \lim_{n\to\infty}f(x_{n})=f(x_{0}) for every f\in X^{*}. The notation for this is x_{n}\stackrel{w}{\to}x_{0}.

Strong convergence(in normed vector space)

||x_{n}-x_{0}||\to0 (n\to\infty)

The notation for this is x_{n}\to x_{0} or x_{n}\stackrel{s}{\to}x_{0}. We also call it convergence in norm.

Weak* convergence

Suppose X is a normed vector space, \{f_{n}\}\subset X^{*}, f_{0}\subset X^{*}(n=1,2,\cdots).

It is the same as strong convergence for a series of operators.

For any x\in X, f_{n}(x)\to f_{0}(x). The notation for this is f_{n}\stackrel{w*}{\to}f_{0}.

Strong convergence

If \{f_{n}\} converse to f_{0} by the norm on X^{*}, or ||f_{n}-f_{0}||\to0(n\to\infty)

It is the same as uniformly convergence for a series of operators.

References:

1. http://en.wikipedia.org

2. Optimization by Vector Space Methods, Luenberger, David G., 1998

Posted in convergence, space | 1 Comment

Understand the Maximum Likelihood Estimation

Let X_{1},X_{2}\cdots, be iid random vector with density or probability function f(x,\theta), where \theta is the unknown parameter and suppose the true value of \theta is \theta_{0}.

The likelihood function: L(x,\theta)=\prod_{i=1}^{n}f(X_{i},\theta)
The logarithm likelihood function: l(x,\theta)=\Sigma_{i=1}^{n}logf(X_{i},\theta)

Note that these two functions are actually joint probability functions of the iid data X_{1},X_{2},\cdots; we call them “likelihood” functions instead of “probability” functions because we now consider them as functions of the unknown parameter \theta.

The maximum likelihood estimation \hat{\theta}_{MLE} is found by maximize the likelihood function, and we will show the idea behind this procedure.

Jessen inequality states that for a concave function g, Eg(X)\leq g(EX) for any random variable X. g(x)=logx is concave, so under \theta_{0}, E_{\theta_{0}}logf(X,\theta)-E_{\theta_{0}}logf(X,\theta_{0})=E_{\theta_{0}}log\frac{f(X,\theta)}{f(X,\theta_{0})}\leq log[E_{\theta_{0}}\frac{f(X,\theta)}{f(X,\theta_{0})}]=log[\int\frac{f(X,\theta)}{f(X,\theta_{0})}f(x,\theta_{0})dx]=log1=0.

That is to say, max_{\theta}E_{\theta_{0}}logf(X,\theta)=E_{\theta_{0}}logf(X,\theta_{0}), since E_{\theta_{0}}logf(X,\theta)\leq E_{\theta_{0}}logf(X,\theta_{0}).

Now a natural criteria for finding \theta_{0} is to find the parameter which maximizes h(\theta)\equiv E_{\theta_{0}}logf(X,\theta).

However, the function we want to maximize h(\theta) is unknown because h(\theta) has something to do with the unknown \theta_{0}.

But note that \theta_{0} only appears in the mean operator E_{\theta_{0}}, we can overcome this problem by using sample mean operator \frac{1}{n}\Sigma instead.

Finally, we reach that we can find a estimator by maximize l(x,\theta)=\Sigma_{i=1}^{n}logf(X_{i},\theta).

There will be a loss when using sample mean to replace the real one, and the quality of this estimator depends on the law of large numbers.

(Update: 2012/Feb/17) We want to maximize logP(X;\theta), or equivalently, letting

(logP(X;\theta))'=\frac{P'(X;\theta)}{P(X;\theta)}\to0.

The above formula can be interpreted as we want

1. P(X;\theta)\to1, or as large as possible,

2. P'(X;\theta)\to0, or as small as possible.

These two points implies that we hope our estimation \hat{\theta} can make the possibility of the observed data X as large as possible, and at the same time, we want some kind of stability of that possibility, i.e., we hope the rate of change in P(X;\theta) – that is P'(X;\theta), as small as possible.

References:
A Course in Large Sample Theory(Lecture notes), Xianyi Wu
A Course in Large Sample Theory, Thomas S. Ferguson

Posted in statistics | Leave a comment

Convergence in Probability and Convergence with Probability 1

Definations:

For a random vector X and a sequence of random vectors X_{1},\cdots X_{n},\cdots defined on a probability space (\Omega,\mathcal{F},P).

1. Convergence in Probability or Convergence of probability measures, is kind of convergence in measure, in concept of the measure theory.
X_{n}\stackrel{P}{\to}X: P\{||X_{n}-X||>\epsilon\}\to0 for any fixed \epsilon>0.

2. Convergence with Probability 1, is kind of pointwise convergence in real analysis, it is also called almost surely convergence in the measure theory.
X_{n}\stackrel{a.s.}{\to}X: if P\{X_{n}\to X\}=1 or P\{X_{n}\nrightarrow X\}=0.

Convergence in Probability is weaker than Convergence with Probability 1, as stated below:

Theorem: X_{n}\stackrel{a.s.}{\to}X\Rightarrow X_{n}\stackrel{P}{\to}X

Proof: \{X_{n}\to X\}={\omega: for each \epsilon>0 there exists an N(\omega)>0 such that ||X_{k}(\omega)-X(\omega)||\leq\epsilon for all k>N(\omega)}

={\cap_{\epsilon>0}{\omega: there exists an N(\omega)>0 such that ||X_{k}(\omega)-X(\omega)||\leq\epsilon for all k>N(\omega)}}

={\cap_{\epsilon>0}\cup_{n=1}^{\infty}{\omega: ||X_{k}(\omega)-X(\omega)||\leq\epsilon for all k>N(\omega)}}

={\cap_{\epsilon>0}\cup_{n=1}^{\infty}\cap_{k>n}{\omega: ||X_{k}(\omega)-X(\omega)||\leq\epsilon}}

So X_{n}\stackrel{a.s.}{\to}X\Longleftrightarrow P\{X_{n}\to X\}=1

\Longleftrightarrow P\{\cap_{\epsilon>0}\cup_{n=1}^{\infty}\cap_{k>n}\{\omega:||X_{k}(\omega)-X(\omega)||\leq\epsilon\}\}=1

\Longleftrightarrow P\{\cup_{n=1}^{\infty}\cap_{k>n}\{\omega:||X_{k}(\omega)-X(\omega)||\leq\epsilon\}\}=1 for all \epsilon>0

\Longleftrightarrow\lim_{n\to\infty}P\{\cap_{k>n}\{\omega:||X_{k}(\omega)-X(\omega)||\leq\epsilon\}\}=1 for all \epsilon>0

\Longleftrightarrow\lim_{n\to\infty}P\{\cup_{k>n}\{\omega:||X_{k}(\omega)-X(\omega)||>\epsilon\}\}=0 for all \epsilon>0

\Rightarrow\lim_{n\to\infty}P\{||X_{n}-X||>\epsilon\}=0 for all \epsilon>0, or P\{||X_{n}-X||>\epsilon\}\to0 for any fixed \epsilon>0, that is X_{n}\stackrel{P}{\to}X.

if X_{n}\stackrel{a.s.}{\to}X, let the set of unconvergence points as U_{a.s.}=\{\omega:X_{n}(\omega)\nrightarrow X(\omega)\} then P\{U_{a.s.}\}=0.
if X_{n}\stackrel{P}{\to}X, let the set of unconvergence points as U_{P}=\{\omega:X_{n}(\omega)\nrightarrow X(\omega)\} then we donot have P\{U_{P}\}=0.

Convergence in Probability do not require P\{U_{P}\}=0 but release the condition to \lim_{n\to\infty}P\{U_{P}(n)\}=0, where U_{P}(n) is a set which will go smaller as n\to\infty.

From the proof above we can see that U_{P}=\{\cup_{\epsilon}\cap_{n=1}^{\infty}\cup_{k}\{\omega:||X_{k}(\omega)-X(\omega)||>\epsilon\}\} and U_{P}(n)=\{\omega:||X_{n}-X||>\epsilon\}, that is to say U_{P} is bigger than U_{a.s.}

References:
A Course in Large Sample Theory(Lecture notes), Xianyi Wu
A Course in Large Sample Theory, Thomas S. Ferguson

Posted in analysis, convergence, probability | 24 Comments

Stochastic o and O symbols

Stochastic o and O symbols are the basic symbols for Asymptotic Statistics or Large Sample Theory.

(i) A_{n}=o_{p}(B_{n}): if |\frac{A_{n}}{b_{n}}|\stackrel{P}{\to}0.

sequence of random variables A_{n} is of smaller order in probability than a sequence B_{n}.

In particular, A_{n}=o_{p}(1), “small oh-P-one”, if and only if A_{n}\stackrel{P}{\to}0; so A_{n}=o_{p}(B_{n}) means A_{n}=Y_{n}B_{n} and Y_{n}\stackrel{P}{\to}0.

Example: X=o_{p}(1) means X\stackrel{P}{\to}0, and X=o_{p}(n^{-1/2}) means n^{1/2}X\stackrel{P}{\to}0, or X goes to 0 faster than \frac{1}{n^{1/2}} in probability(such as X=\frac{1}{n}).

(ii) A_{n}=O_{p}(B_{n}) : if given \epsilon>0, there exists a constant M=M(\epsilon) and an integer n_{0}=n_{0}(\epsilon) such that P(|A_{n}|\leq M|B_{n}|)\geq1-\epsilon for all n>n_{0}.

sequence A_{n} is be of order less than or equal to that of B_{n} in probability.

In particular, A_{n}=O{}_{p}(1), “big oh-P-one”, if for any\epsilon>0, there exists a constant M and an integer n_{0} such that P(|A_{n}|\leq M)\geq1-\epsilon for all n>n_{0}, A_{n} is said to be bounded in probability(or tight); so A_{n}=O_{p}(B_{n}) means A_{n}=Y_{n}B_{n} and Y_{n}=O_{p}(1).

It’s easy to see from the definition that O_{p}(1)=O_{p}(C) for any constant 0<C<\infty.

(iii) A_{n}\asymp_{p}B_{n}: if given \epsilon>0, there exist constants 0<m<M<\infty and an integer n_{0} such that P[m<|\frac{A_{n}}{B_{n}}|<M]\geq1-\epsilon for all n>n_{0}.

sequence A_{n} is said to be of the same order as B_{n} in probability.

Some facts:

o_{p}(1)+o_{p}(1)=o_{p}(1): If X_{n}\stackrel{P}{\to}0 and Y_{n}\stackrel{P}{\to}0, then Z_{n}=X_{n}+Y_{n}\stackrel{P}{\to}0. (example of continuous-mapping theorem)

o_{p}(1)+O_{p}(1)=O_{p}(1)

O_{p}(1)o_{p}(1)=o_{p}(1): If the sequence \{Y_{n},n=1,2,\cdots\} is bounded in probability and if \{C_{n}\} is a sequence of random variables tending to 0 in probability, then C_{n}Y_{n}\stackrel{P}{\to}0.

(1+o_{p}(1))^{-1}=O_{p}(1)

o_{p}(R_{n})=R_{n}o_{p}(1)

O_{p}(R_{n})=R_{n}O_{p}(1)

o_{p}(O_{p}(1))=o_{p}(1)

Lemma: Let R be a function defined on domain in \mathcal{R}^{k} such that R(0)=0. Let X_{n} be a sequence of random vectors with values in the domain of R that converges in probability to zero. Then, for every p>0,

(i) if R(h)=o(||h||^{p}) as h\to0, then R(X_{n})=o_{p}(||X_{n}||^{p});

(ii) if R(h)=O(||h||^{p}) as h\to0, then R(X_{n})=O_{p}(||X_{n}||^{p});

Result: For a random variable S, S=ES+O_{p}(\sqrt{Var(S})).

Proof:

We only needs to prove that (S-ES)/\sqrt{Var(S)}=O_{p}(1) or equally, for any\epsilon>0, there exists a constant M and an integer n_{0} such that P(|(S-ES)/\sqrt{Var(S)}|\leq M)\geq1-\epsilon for all n>n_{0}.

Let NS=(S-ES)/\sqrt{Var(S)},

According to Markov inequality, P(|NS|\leq M)\geq ENS^{2}/M^{2}=[E(S-ES)^{2}/Var(S)]/M^{2}=1/M^{2}\to0 as M\to\infty.

From the proof above we know that for any normalized random variable NS=(S-ES)/\sqrt{Var(S)}, we have NS=O_{p}(1), or NS is bounded in probability- the reason is natural, if any random variable is not bounded, either its mean is too large(E(S_{n})\to\infty) or it varies too much(Var(S_{n})\to\infty), and normalization will eliminate those two possibilities. On the other hand, for a specified random variable S, if ES<\infty and Var(S)<\infty, then S=ES+O_{p}(1), especially, when ES=0, S=O_{p}(1) .

Example: from center limit theorem we know that \sqrt{n}(\bar{X}-EX)\to N(0,DX), then we have

//\sqrt{n}(\bar{X}-EX)=N(0,DX)+o_{p}(1)=\sqrt{DX}O_{p}(1)+o_{p}(1)=O_{p}(1)//,

//\bar{X}=EX+O_{p}(1)\times n^{-1/2}=EX+O_{p}(n^{-1/2})//.

P(\frac{\sqrt{n}(\bar{X}-EX)}{DX}>M )\to P(Z>M),\ Z\sim N(0,1). P(Z>M) can be smaller than \forall\epsilon as long as M is large enough, so \frac{\sqrt{n}(\bar{X}-EX)}{DX}=O_p(1), or \bar{X}=EX+O_{p}(n^{-1/2}).

The weak law of large numbers states that \bar{X}\stackrel{P}{\to}EX, so we have

\bar{X}-EX=o_{p}(1).

(Update: 2012/Feb/17) Similarly, let X_{n} be a sequence of random vectors, using Markov inequality P(|X_{n}|>M)\leq\frac{E|X_{n}|^{k}}{M^{k}}, we have

1. If there is a number k>0 such that E|X_{n}|^{k} is bounded, then X_{n}=O_{p}(1);
similarly, if E|X_{n}|^{k}\leq ca_{n}, where c is a constant and a_{n} is a sequence of positive numbers,
then X_{n}=O_{p}(a_{n}^{1/k}).

2. If there is a number k>0 such that E|X_{n}|^{k}\to0 (So M can be \epsilon), then X_{n}=o_{p}(1);
similarly, if E|X_{n}|^{k}\leq ca_{n}, where c is a constant and a_{n} is a sequence of positive numbers,
then X_{n}=o_{p}(b_{n}) for any sequence b_{n}>0 such that b_{n}^{-1}a_{n}^{1/k}\to0.

3. If there are sequences of vectors \{\mu_{n}\} and singularization matrices \{A_{n}\} such that A_{n}(X_{n}-\mu_{n}) converges in distribution, then X_{n}=\mu_{n}+O_{p}(||A_{n}^{-1}||).

References:

Elements of Large-Sample Theory, E.L. Lehmann, 1998

Asymptotic Statistics, A. W. van der Vaart, 2000

Linear and Generalized Linear Mixed Models and Their Applications, Jiming Jiang, 2006

Posted in analysis, asymptotic, convergence | Leave a comment

Using Lib FANN in R via Rcpp

Test Environment: OS: Windows XP, R: R 2.12.2(G:\Program Files\R\R-2.12.2\), Rcpp: 0.9.4, Rtools: 212 (G:\Rtools), FANN: 2.1.0
Lib FANN(Fast Artificial Neural Network Library) is a free open source neural network library written in C. In this tutorial I will show you how to use the Rcpp package to wrap FANN to R.

1. Download the FANN library: http://leenissen.dk/fann/wp/download/ (fann-2.1.0beta.zip)

2. Unzip the source code to some folder(“G:\libFANN\fann-2.1.0\“)

3. Make sure that you installed Rtools and setted the path for MinGW properly. You can check this post for that. The path for my MinGW(Rtools) is “G:\Rtools\MinGW

4. Install MSYS: Download and install MSYS 1.0.11 from http://www.mingw.org/wiki/MSYS. I installed MSYS into “G:\msys“. Make sure to tell msys where your installed MinGW(Rtools) is. This can be done in the post-install process of MSYS, or manually create a file named “fstab” in folder “G:\msys\1.0\etc“, with the following line:
G:/Rtools/MinGW /mingw #tell msys the path of MinGW

5. Run msys(there should be a shortcut on your desktop: “G:\msys\1.0\msys.bat -norxvt“), do the config for libFANN:
cd G:\libFANN\fann-2.1.0
./configure
Now you should find the Makefile created in “G:\libFANN\fann-2.1.0\src, let’s compile the source code of FANN:
cd src
make
After the compilation, you will find some files such as “fixedfann.o, floatfann.o, doublefann.o” created in the src folder.

6. Let’s first test the compiled lib: Copy the files “xor.data, xor_train.c, xor_test.c” from “examples” folder to “src” folder.
ar rc libFANNLib.a doublefann.o fixedfann.o floatfann.o
ranlib libFANNLib.a
gcc -o xor_train xor_train.c -Iinclude libFANNLib.a
gcc -o xor_test xor_test.c -Iinclude libFANNLib.a
Now will can find “xor_train.exe, xor_test.exe” in the “src” folder, click one after the other to run the test(to view the results, you should run them both from cmd).

7. It’s time to wrap FANN to R, we will write a simple package “RcppFANN” using Rcpp Module: go to R

library(Rcpp);
Rcpp.package.skeleton("RcppFANN");

find the created files in “C:\Documents and Settings\Administrator\My Documents\RcppFANN“(your working directory)

8. Go to the folder “C:\Documents and Settings\Administrator\My Documents\RcppFANN\src” and edit the two files “rcpp_hello_world.h, rcpp_hello_world.cpp“, we simply wrap the example code “xor_train.c, xor_test.c” into “rcpp_hello_world.cpp” using Rcpp Module:

//rcpp_hello_world.h
#ifndef _RcppFANN_RCPP_HELLO_WORLD_H
#define _RcppFANN_RCPP_HELLO_WORLD_H

#include

/*
 * note : RcppExport is an alias to `extern "C"` defined by Rcpp.
 *
 * It gives C calling convention to the rcpp_hello_world function so that
 * it can be called from .Call in R. Otherwise, the C++ compiler mangles the
 * name of the function and .Call can't find it.
 *
 * It is only useful to use RcppExport when the function is intended to be called
 * by .Call. See the thread http://thread.gmane.org/gmane.comp.lang.r.rcpp/649/focus=672
 * on Rcpp-devel for a misuse of RcppExport
 */

#endif

 

//rcpp_hello_world.cpp
#include "rcpp_hello_world.h"

using namespace Rcpp ;

#include

#include "fann.h"

int xor_test()
{
	fann_type *calc_out;
	unsigned int i;
	int ret = 0;

	struct fann *ann;
	struct fann_train_data *data;

	printf("Creating network.\n");

#ifdef FIXEDFANN
	ann = fann_create_from_file("xor_fixed.net");
#else
	ann = fann_create_from_file("xor_float.net");
#endif

	if(!ann)
	{
		printf("Error creating ann --- ABORTING.\n");
		return 0;
	}

	fann_print_connections(ann);
	fann_print_parameters(ann);

	printf("Testing network.\n");

#ifdef FIXEDFANN
	data = fann_read_train_from_file("xor_fixed.data");
#else
	data = fann_read_train_from_file("xor.data");
#endif

	for(i = 0; i < fann_length_train_data(data); i++) 	{ 		fann_reset_MSE(ann); 		calc_out = fann_test(ann, data->input[i], data->output[i]);
#ifdef FIXEDFANN
		printf("XOR test (%d, %d) -> %d, should be %d, difference=%f\n",
			   data->input[i][0], data->input[i][1], calc_out[0], data->output[i][0],
			   (float) fann_abs(calc_out[0] - data->output[i][0]) / fann_get_multiplier(ann));

		if((float) fann_abs(calc_out[0] - data->output[i][0]) / fann_get_multiplier(ann) > 0.2)
		{
			printf("Test failed\n");
			ret = -1;
		}
#else
		printf("XOR test (%f, %f) -> %f, should be %f, difference=%f\n",
			   data->input[i][0], data->input[i][1], calc_out[0], data->output[i][0],
			   (float) fann_abs(calc_out[0] - data->output[i][0]));
#endif
	}

	printf("Cleaning up.\n");
	fann_destroy_train(data);
	fann_destroy(ann);

	return ret;
}

////

int FANN_API test_callback(struct fann *ann, struct fann_train_data *train,
	unsigned int max_epochs, unsigned int epochs_between_reports,
	float desired_error, unsigned int epochs)
{
	printf("Epochs     %8d. MSE: %.5f. Desired-MSE: %.5f\n", epochs, fann_get_MSE(ann), desired_error);
	return 0;
}

int xor_train()
{
	fann_type *calc_out;
	const unsigned int num_input = 2;
	const unsigned int num_output = 1;
	const unsigned int num_layers = 3;
	const unsigned int num_neurons_hidden = 3;
	const float desired_error = (const float) 0;
	const unsigned int max_epochs = 1000;
	const unsigned int epochs_between_reports = 10;
	struct fann *ann;
	struct fann_train_data *data;

	unsigned int i = 0;
	unsigned int decimal_point;

	printf("Creating network.\n");
	ann = fann_create_standard(num_layers, num_input, num_neurons_hidden, num_output);

	data = fann_read_train_from_file("xor.data");

	fann_set_activation_steepness_hidden(ann, 1);
	fann_set_activation_steepness_output(ann, 1);

	fann_set_activation_function_hidden(ann, FANN_SIGMOID_SYMMETRIC);
	fann_set_activation_function_output(ann, FANN_SIGMOID_SYMMETRIC);

	fann_set_train_stop_function(ann, FANN_STOPFUNC_BIT);
	fann_set_bit_fail_limit(ann, 0.01f);

	fann_init_weights(ann, data);

	printf("Training network.\n");
	fann_train_on_data(ann, data, max_epochs, epochs_between_reports, desired_error);

	printf("Testing network. %f\n", fann_test_data(ann, data));

	for(i = 0; i < fann_length_train_data(data); i++) 	{ 		calc_out = fann_run(ann, data->input[i]);
		printf("XOR test (%f,%f) -> %f, should be %f, difference=%f\n",
			   data->input[i][0], data->input[i][1], calc_out[0], data->output[i][0],
			   fann_abs(calc_out[0] - data->output[i][0]));
	}

	printf("Saving network.\n");

	fann_save(ann, "xor_float.net");

	decimal_point = fann_save_to_fixed(ann, "xor_fixed.net");
	fann_save_train_to_fixed(data, "xor_fixed.data", decimal_point);

	printf("Cleaning up.\n");
	fann_destroy_train(data);
	fann_destroy(ann);

	return 0;
}

////

RCPP_MODULE(mod){
 function("xor_test", &xor_test);
 function("xor_train", &xor_train);
}

9. Edit the “Makevars.win” to:
## Use the R_HOME indirection to support installations of multiple R version
PKG_LIBS = $(shell “${R_HOME}/bin${R_ARCH_BIN}/Rscript.exe” -e “Rcpp:::LdFlags()”) G:/libFANN/fann-2.1.0/src/libFANNLib.a
PKG_CPPFLAGS = -I”G:\libFANN\fann-2.1.0\src\include”

10. Now we can compile and install our test package:
Windows Start -> run -> cmd -> cd C:\Documents and Settings\Administrator\My Documents
C:\Documents and Settings\Administrator\My Documents>
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” check RcppFANN
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” INSTALL RcppFANN

11. Finally, let test our package:
copy the data files “xor.data“(located in “G:\libFANN\fann-2.1.0\examples“) to folder “C:\Documents and Settings\Administrator\My Documents\“(or your working directory).
In R:

library("Rcpp");
library("RcppFANN");
mod mod$xor_train();

Result:

[1] 0

If everything goes well, you can find three files “xor_fixed.data, xor_fixed.net, xor_float.net“(the trained network) created in “C:\Documents and Settings\Administrator\My Documents\“.

mod$xor_test();

Result:

[1] 0

Well, FANN is successfully running in R, but we can’t see anything useful in R because the results of the original test function are printed using the C “printf” function. To make the lib really useful, you can write some further wrapper functions to retrieve those results and supply data to FANN in R.

Source code of the test package “RcppFANN”: https://bzstat.googlecode.com/svn/trunk/blogpost/RcppFANN

Links:
http://dirk.eddelbuettel.com/code/rcpp.html
http://dirk.eddelbuettel.com/code/rcpp/Rcpp-package.pdf
http://dirk.eddelbuettel.com/code/rcpp/Rcpp-modules.pdf
http://stat.ethz.ch/R-manual/R-devel/doc/manual/R-admin.html
http://sourceforge.net/projects/rfann/
http://cran.r-project.org/web/packages/RSNNS/

Posted in programming, r | Leave a comment

The Rcpp package test notes

Test Environment: OS: Windows XP, R: R 2.12.2(G:\Program Files\R\R-2.12.2\), Rcpp: 0.9.4

1. Download and install the Rtools(for windows): http://www.murdoch-sutherland.com/Rtools/
I installed the Rtools into “G:\Rtools“.

2. Set Enviroment Path for windows:
Desktop -> right-click My Computer -> click Properties -> click the Advanced tab -> click Environment Variables button
-> highlight the “Path” variable in the Systems Variable section -> click the Edit button -> addG:\Rtools\bin;G:\Rtools\MinGW\bin;” there (“C:\MiKTeX\miktex\bin;” may needed too.)
and
in the Systems Variable section -> click the New button -> add a new variable named “CYGWIN” with value “nodosfilewarning” (Turn off MS-DOS style path warnings with cygwin utilities)

3. Install the R package “Rcpp” and “inline”.

4. Quick test using inline

library(inline);
library("Rcpp");

inc <- '
double norm(double x, double y){
return sqrt(x*x+y*y);
}

RCPP_MODULE(mod){
 function("norm", &norm);
}
'
fx <- cxxfunction( signature(), "" , include = inc, plugin = "Rcpp" );
mod <- Module( "mod", getDynLib(fx) );
mod$norm(3,4);

Result:

[1] 5

5. Quick test of writing a simple package

Rcpp.package.skeleton("RcppTestPackage");
Creating directories ...
Creating DESCRIPTION ...
Creating NAMESPACE ...
Creating Read-and-delete-me ...
Saving functions and data ...
Making help files ...
Done.
Further steps are described in './RcppTestPackage/Read-and-delete-me'.

Adding Rcpp settings
 >> added Depends: Rcpp
 >> added LinkingTo: Rcpp
 >> added useDynLib directive to NAMESPACE
 >> added Makevars file with Rcpp settings
 >> added Makevars.win file with Rcpp settings
 >> added example header file using Rcpp classes
 >> added example src file using Rcpp classes
 >> added example R file calling the C++ example
 >> added Rd file for rcpp_hello_world

Now you can find the created file in “C:\Documents and Settings\Administrator\My Documents\RcppTestPackage” (or your working directory)

Now we can compile and install the test package:
Windows Start -> run -> cmd -> cd C:\Documents and Settings\Administrator\My Documents
C:\Documents and Settings\Administrator\My Documents>
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” check RcppTestPackage
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” INSTALL RcppTestPackage
(“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” build RcppTestPackage (This is the file you can upload to CRAN.))

Finally, go to R and test your package:

library(RcppTestPackage);
rcpp_hello_world();

Result:

[[1]]
[1] "foo" "bar"

[[2]]
[1] 0 1

6. Create a simple package using Rcpp Module:

Rcpp.package.skeleton("RcppTestModule");

Find the created file in “C:\Documents and Settings\Administrator\My Documents\RcppTestPacakge” (or your working directory)
Go to the “src” folder and make changes to the files “rcpp_hello_world.h” and “rcpp_hello_world.cpp”,
Let’s write something meaningful:

//rcpp_hello_world.h
#ifndef _RcppTestPacakge_RCPP_HELLO_WORLD_H
#define _RcppTestPacakge_RCPP_HELLO_WORLD_H

#include <Rcpp.h>

/*
 * note : RcppExport is an alias to `extern "C"` defined by Rcpp.
 *
 * It gives C calling convention to the rcpp_hello_world function so that 
 * it can be called from .Call in R. Otherwise, the C++ compiler mangles the 
 * name of the function and .Call can't find it.
 *
 * It is only useful to use RcppExport when the function is intended to be called
 * by .Call. See the thread http://thread.gmane.org/gmane.comp.lang.r.rcpp/649/focus=672
 * on Rcpp-devel for a misuse of RcppExport
 */

#endif
//rcpp_hello_world.cpp
#include "rcpp_hello_world.h"

using namespace Rcpp;

double norm(double x, double y){
return sqrt(x*x+y*y);
}

RCPP_MODULE(mod){
 function("norm", 
	      &norm,
		  List::create( _["x"] = 0.0, _["y"] = 0.0),
		  "Provides a simple vector norm");
}

Similarly, compile and install the test package:
Windows Start -> run -> cmd -> cd C:\Documents and Settings\Administrator\My Documents
C:\Documents and Settings\Administrator\My Documents>
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” check RcppTestModule
“G:\Program Files\R\R-2.12.2\bin\i386\Rcmd” INSTALL RcppTestModule

Finally, test our package in R:

library("RcppTestModule");
mod <- Module( "mod", PACKAGE = "RcppTestModule");
mod$norm(3,4);

Result:

[1] 5

and

mod$norm();

Result:

[1] 0

also

show(mod$norm);

Result:

internal C++ function <0241A0A8>
    docstring : Provides a simple vector norm
    signature : double norm(double, double)

Source code of package “RcppTestModule”: https://bzstat.googlecode.com/svn/trunk/blogpost/RcppTestModule

Links:

http://dirk.eddelbuettel.com/code/rcpp.html
http://dirk.eddelbuettel.com/code/rcpp/Rcpp-package.pdf
http://dirk.eddelbuettel.com/code/rcpp/Rcpp-modules.pdf
http://cran.r-project.org/doc/manuals/R-exts.html
http://www.stat.psu.edu/~dsy109/SOS_Talk.pdf
http://cran.r-project.org/doc/contrib/Leisch-CreatingPackages.pdf

Posted in programming, r | 1 Comment

R Tricks – Macros and Conditional Execution in R

1 Use Macros in R

Problem

You want to use C-like macros, such as “#define” in R.

Solution

Use “defmacro” in {gtools} pacakge.

> library(gtools);
> mul<-defmacro(a, b, expr={a*b});
> mul(2,3);
[1] 6
> CONST_PI<-defmacro(0, expr={3.1415926});#0 is useless here
> print(CONST_PI(0));
[1] 3.141593

http://cran.r-project.org/web/packages/gtools/index.html
https://stat.ethz.ch/pipermail/r-help/2007-April/130015.html

2 Conditional Execution

Problem

You want to use C-like macros, such as “#include” and C-like conditional compiling macros, such as “#if, #elif, #else, #endif“, in R.

Solution

There is no preprocessor for R, but you can use keyword “source” and “if/else” for conditional execution.

Debug<-TRUE;
if(Debug)
{source("C:/Debug.r");}
else
{source("C:/NoDebug.r");}
#C:\Debug.r
print("Debug");
#C:\NoDebug.r
print("NoDebug");
Posted in programming, r | Leave a comment