diff --git a/R/pensionTable.R b/R/pensionTable.R
index 5511b1daa5e9a6cb23abf3cdc97981efc1a802f9..7b8913aa2467be6566662c748d6f840d404eeaca 100644
--- a/R/pensionTable.R
+++ b/R/pensionTable.R
@@ -28,7 +28,7 @@ NULL
 #'
 #' @slot qx     Death probability table of actives (derived from mortalityTable)
 #' @slot ix     Invalidity probability of actives (derived from mortalityTable)
-#' @slot qxi    Death probability table of invalids (derived from mortalityTable)
+#' @slot qix    Death probability table of invalids (derived from mortalityTable)
 #' @slot rx     Reactivation probability of invalids (derived from mortalityTable)
 #' @slot apx    Retirement probability of actives (derived from mortalityTable)
 #' @slot apix   Retirement probability of invalids (derived from mortalityTable)
@@ -66,17 +66,18 @@ pensionTable = setClass(
 #'
 #' @examples
 #' pensionTables.load("Austria_*", wildcard=TRUE)
-#' transitionProbabilities(EttlPagler.male)
+#' # transitionProbabilities(EttlPagler.male)
 #'
 #' @exportMethod transitionProbabilities
-setGeneric("transitionProbabilities", function(object, ...) standardGeneric("transitionProbabilities"));
+setGeneric("transitionProbabilities", function(object, ..., YOB = 1982) standardGeneric("transitionProbabilities"));
 
-#' @describeIn baseTable Return the base table of the joint lives mortality table (returns the base table of the first table used for joint lives)
+#' @describeIn transitionProbabilities Return all transition probabilities of the pension table
 setMethod("transitionProbabilities", "pensionTable",
-          function(object, ...,  YOB = 1982) {
+          function(object, ..., as.data.frame = TRUE,  YOB = 1982) {
+              na.zero = function(x) { x[is.na(x)] = 0; x }
               x = ages(object@qx);
-              q = deathProbabilities(object@qx, ..., YOB = YOB);
-              i = deathProbabilities(object@ix, ..., YOB = YOB);
+              q = na.zero(deathProbabilities(object@qx, ..., YOB = YOB));
+              i = na.zero(deathProbabilities(object@ix, ..., YOB = YOB));
               qi = deathProbabilities(object@qix, ..., YOB = YOB);
               r = deathProbabilities(object@rx, ..., YOB = YOB);
               ap = deathProbabilities(object@apx, ..., YOB = YOB);
@@ -86,44 +87,144 @@ setMethod("transitionProbabilities", "pensionTable",
               qw = deathProbabilities(object@qwy, ..., YOB = YOB);
               yx = deathProbabilities(object@yx, ..., YOB = YOB);
               qg = deathProbabilities(object@qgx, ..., YOB = YOB);
-              data.frame(x, q, i, qi, r, ap, api, qp, h, qw, yx, qg)
+              if (as.data.frame) {
+                  data.frame(x, q, i, qi, r, ap, api, qp, h, qw, yx, qg)
+              } else {
+                  states = c("a", "i", "p", "d")
+                  transProb = array(0, dim = c(4,4, length(x)), dimnames = list(states, states, x))
+
+                  transProb["a", "a", ] = (1 - i - q) * (1 - ap);
+                  transProb["a", "i", ] = i;
+                  transProb["a", "p", ] = (1 - q - i ) * ap;
+                  transProb["a", "d", ] = q;
+
+                  transProb["i", "a", ] = r;
+                  transProb["i", "i", ] = (1 - qi - r) * (1 - api);
+                  transProb["i", "p", ] = (1 - qi - r) * api;
+                  transProb["i", "d", ] = qi;
+
+                  transProb["p", "p", ] = 1 - qp;
+                  transProb["p", "d", ] = qp;
+
+                  transProb["d", "d", ] = 1;
+
+                  list(transitionProbabilities = transProb, widows = data.frame(x, h, qw, yx))
+              }
           })
 
 
 if (FALSE) {
+    transitionProbabilities(AVOe2008P.male, YOB = 1977, as.data.frame = FALSE)
     epP = transitionProbabilities(EttlPagler.male, YOB = 1982)
-    avoe08p = transitionProbabilities(AVOe2008P.male, YOB = 1977)
+#    avoe08p =
+        transitionProbabilities(AVOe2008P.male, YOB = 1977, as.data.frame = TRUE)
+}
+
+bwRente = function(p, v) {
+    Reduce(function(pp, ax1) { 1 + pp * ax1 * v }, p, 0.0, right = TRUE, accumulate = TRUE)[-(length(p) + 1)];
 }
 
 
+reservesThieleRecursion = function(p, ai, aij, states, i = 0.03) {
+    v = 1 / (1 + i)
+    res = array(0, dim = dim(ai), dimnames = dimnames(ai));
+    # Recursive relation:
+    #    Vi(t,A) = ai(t) + \sum_j v p_ij(t) (aij(t) + Vj(t+1,A))
+    # with: ai(t) .. payment at t for being in state i
+    #       aij(t) ... payment at t+1 for switching from state i to j
+    #       Vi(t,A) ... reserve for payments A in state i at time t
+    ThieleRecursion = function(t, Vt1) {
+        rr = ai[,t] + v * rowSums(p[,,t] * aij[,,t]) + v * as.vector(p[,,t] %*% Vt1)
+        as.vector(rr)
+    }
+    # Loop backwards over all times (starting value for reserves is 0)
+    times = dimnames(p)[[3]];
+    res = Reduce(f = ThieleRecursion, x = times, init = rep(0, length(states)), right = TRUE, accumulate = TRUE)[-(length(times) + 1)]
+    res = do.call("cbind", res)
+    dimnames(res) = dimnames(ai)
+    res
+}
+if (FALSE) {
+res = anwartschaften(AVOe2008P.female, YOB = 1977);
+res
+}
+
+
+#' Calculates all "anwartschaften" for the gien pension table
+#'
+#' @param object A pension table object (instance of a \code{\linkS4class{pensionTable}} class)
+#' @param ... Currently unused
+#' @param i Interest rate (default 0.03)
+#' @param YOB Year of birth (default 1982)
+#'
+#' @examples
+#' pensionTables.load("Austria_*", wildcard=TRUE)
+#' # anwartschaften(EttlPagler.male, i=0.03, YOB=1972)
+#'
+#' @exportMethod transitionProbabilities
 setGeneric("anwartschaften", function(object, ...) standardGeneric("anwartschaften"));
 
+#' @describeIn anwartschaften Calculates all "anwartschaften" for the gien pension table
 setMethod("anwartschaften", "pensionTable",
-    function(object, ...,  i = 0.03, YOB = 1982) {
-        probs = transitionProbabilities(object, ..., YOB);
-        anwartschaften(probs, ..., YOB)
-    }
-);
+          function(object, ...,  i = 0.03, YOB = 1982) {
+              probs = transitionProbabilities(object, ..., YOB = YOB, as.data.frame = FALSE);
 
-bwRente = function(p, v) {
-    Reduce(function(pp, ax1) { 1 + pp * ax1 * v }, p, 0.0, right = TRUE, accumulate = TRUE)[-(length(p) + 1)];
-}
+              # Time series of transition probabilities
+              pp = probs$transitionProbabilities;
+              x = dimnames(pp)[[3]]
 
-setMethod("anwartschaften", "data.frame",
-    function(object, ...,  i = 0.03) {
-        x = object$x;
-        v = 1 / (1 + i);
-        # Anwartschaft auf Witwenrente und Alterspension
-        # 1) Barwerte:
-        aa = bwRente(1.0 - object$q, v);
-        ai = bwRente(1. - object$qi - object$r, v);
-        ap = bwRente(1. - object$qp, v);
-        aw = bwRente(1. - object$qw, v);
-        data.frame(x, aa, ai, ap, aw)
-    }
-)
+              # Use a data.frame for the annuity PV with the actual ages as dimnames,
+              aw = data.frame(aw = bwRente(1 - probs$widows["qw"], 1 / (1 + i)));
+              dimnames(aw)[[1]] = x
+
+              # Expected death benefit (widows)
+              # Use avg. age of widow to extract the corresponding annuity present value
+              # We used the age as dimname, so we can use simple subsetting
+              expDeathBenefit = probs$widows[["h"]] * aw[as.character(probs$widows[["yx"]]),]
+
+              # Build the matrix of transition payments (only on death there is
+              # the widow PV as benefit, all other transitions do not yield any benefit)
+              states = c("a", "i", "p", "d")
+              transPayments = array(0, dim = c(4,4, length(x)), dimnames = list(states, states, x))
+              transPayments["a","d",] = expDeathBenefit;
+              transPayments["i","d",] = expDeathBenefit;
+              transPayments["p","d",] = expDeathBenefit;
+
+              statePayments = array(0, dim = c(4, length(x)), dimnames = list(states, x));
 
+              aPay = reservesThieleRecursion(p = pp, ai = statePayments + c(1,0,0,0), aij = transPayments*0, states = states, i = i)
+              iPay = reservesThieleRecursion(p = pp, ai = statePayments + c(0,1,0,0), aij = transPayments*0, states = states, i = i)
+              pPay = reservesThieleRecursion(p = pp, ai = statePayments + c(0,0,1,0), aij = transPayments*0, states = states, i = i)
+              wPay = reservesThieleRecursion(p = pp, ai = statePayments, aij = transPayments, states = states)
+
+              list(pp = pp, transPayments = transPayments, statePayments = statePayments, aPay = aPay, iPay = iPay, pPay = pPay, wPay = wPay)
+              list("a" = aPay, "i" = iPay, "p" = pPay, "w" = wPay)
+});
 if (FALSE) {
-    probs = transitionProbabilities(AVOe2008P.female, YOB = 1977)
-    an = anwartschaften(probs, YOB = 1977); an
+res = anwartschaften(AVOe2008P.female, YOB = 1977);
+res
+
+as.array(res$aPay)
+str(res$aPay)
+
+
+dimnames(res$pp)[[3]]
+
+    res["102",]
+res[,"aw"]
+    a=15:43
+a
+a=array(1:8, dim=c(2,4), dimnames=list(c("a1", "a2"), c("b1", "b2", "b3", "b4"))); a
+b=array(11:18, dim=c(2,4), dimnames=list(c("a1", "a2"), c("b1", "b2", "b3", "b4"))); b
+
+array(a, b)
+dimnames(a) = c(15:43)
+
+an = anwartschaften(probs, YOB = 1977); an
+showMethods("anwartschaften")
+showMethods("transitionProbabilities")
+
+
+array(1:12, dim = c(2,3,4), dimnames=list(c("a1", "a2"), c("b1", "b2", "b3"), c("c1", "c2", "c3", "c4")))
+
 }
diff --git a/man/anwartschaften.Rd b/man/anwartschaften.Rd
new file mode 100644
index 0000000000000000000000000000000000000000..6aad13e7b5a4745b90ded17c59ee5dfe519a3846
--- /dev/null
+++ b/man/anwartschaften.Rd
@@ -0,0 +1,34 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/pensionTable.R
+\docType{methods}
+\name{anwartschaften}
+\alias{anwartschaften}
+\alias{anwartschaften,pensionTable-method}
+\title{Calculates all "anwartschaften" for the gien pension table}
+\usage{
+anwartschaften(object, ...)
+
+\S4method{anwartschaften}{pensionTable}(object, ..., i = 0.03, YOB = 1982)
+}
+\arguments{
+\item{object}{A pension table object (instance of a \code{\linkS4class{pensionTable}} class)}
+
+\item{...}{Currently unused}
+
+\item{i}{Interest rate (default 0.03)}
+
+\item{YOB}{Year of birth (default 1982)}
+}
+\description{
+Calculates all "anwartschaften" for the gien pension table
+}
+\section{Methods (by class)}{
+\itemize{
+\item \code{pensionTable}: Calculates all "anwartschaften" for the gien pension table
+}}
+
+\examples{
+pensionTables.load("Austria_*", wildcard=TRUE)
+# anwartschaften(EttlPagler.male, i=0.03, YOB=1972)
+
+}
diff --git a/man/baseTable.Rd b/man/baseTable.Rd
index 01f0f7fe0071f97338b36b3b35efbde35282134b..98c35374d17b90aa96b5dbf937552d51190e593d 100644
--- a/man/baseTable.Rd
+++ b/man/baseTable.Rd
@@ -1,13 +1,11 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/baseTable.R, R/mortalityTable.jointLives.R,
-%   R/pensionTable.R
+% Please edit documentation in R/baseTable.R, R/mortalityTable.jointLives.R
 \docType{methods}
 \name{baseTable}
 \alias{baseTable}
 \alias{baseTable,mortalityTable-method}
 \alias{baseTable,mortalityTable.period-method}
 \alias{baseTable,mortalityTable.jointLives-method}
-\alias{transitionProbabilities,pensionTable-method}
 \title{Return the base table of the life table}
 \usage{
 baseTable(object, ...)
@@ -17,8 +15,6 @@ baseTable(object, ...)
 \S4method{baseTable}{mortalityTable.period}(object, ...)
 
 \S4method{baseTable}{mortalityTable.jointLives}(object, ...)
-
-\S4method{transitionProbabilities}{pensionTable}(object, ..., YOB = 1982)
 }
 \arguments{
 \item{object}{The life table object (class inherited from mortalityTable)}
@@ -35,8 +31,6 @@ Return the base table of the life table
 \item \code{mortalityTable.period}: Return the base table of the life table
 
 \item \code{mortalityTable.jointLives}: Return the base table of the joint lives mortality table (returns the base table of the first table used for joint lives)
-
-\item \code{pensionTable}: Return the base table of the joint lives mortality table (returns the base table of the first table used for joint lives)
 }}
 
 \examples{
diff --git a/man/pensionTable-class.Rd b/man/pensionTable-class.Rd
index 6b8e0989d105586f64753b4207f1077a2e5ef061..d5e8a8312b7d9f1802a9d8a5ef2177a319ca3090 100644
--- a/man/pensionTable-class.Rd
+++ b/man/pensionTable-class.Rd
@@ -36,7 +36,7 @@ Correspondingly, the following transition probabilities can be given:
 
 \item{\code{ix}}{Invalidity probability of actives (derived from mortalityTable)}
 
-\item{\code{qxi}}{Death probability table of invalids (derived from mortalityTable)}
+\item{\code{qix}}{Death probability table of invalids (derived from mortalityTable)}
 
 \item{\code{rx}}{Reactivation probability of invalids (derived from mortalityTable)}
 
diff --git a/man/transitionProbabilities.Rd b/man/transitionProbabilities.Rd
index 843d3823b77076c65596fb7fd3b3c98d33d99576..0b079bccd28e2fce1a36446e0eb968b4b157dd96 100644
--- a/man/transitionProbabilities.Rd
+++ b/man/transitionProbabilities.Rd
@@ -1,10 +1,15 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/pensionTable.R
+\docType{methods}
 \name{transitionProbabilities}
 \alias{transitionProbabilities}
+\alias{transitionProbabilities,pensionTable-method}
 \title{Return all transition probabilities of the pension table}
 \usage{
-transitionProbabilities(object, ...)
+transitionProbabilities(object, ..., YOB = 1982)
+
+\S4method{transitionProbabilities}{pensionTable}(object, ...,
+  as.data.frame = TRUE, YOB = 1982)
 }
 \arguments{
 \item{object}{A pension table object (instance of a \code{\linkS4class{pensionTable}} class)}
@@ -16,8 +21,13 @@ transitionProbabilities(object, ...)
 \description{
 Return all transition probabilities of the pension table
 }
+\section{Methods (by class)}{
+\itemize{
+\item \code{pensionTable}: Return all transition probabilities of the pension table
+}}
+
 \examples{
 pensionTables.load("Austria_*", wildcard=TRUE)
-transitionProbabilities(EttlPagler.male)
+# transitionProbabilities(EttlPagler.male)
 
 }