Package 'tdaR'

Description

Marks a set of columns, used inside a tda_qreg formula, as the same underlying variable's value at each category of the response – price for choice 1, choice 2, choice 3, say, one coefficient shared across all of them – rather than the formula's ordinary one-column-per-predictor, shared-across-categories shape. TDA's "(Z1,Z2,...)" grouped variable-list syntax (what its source calls "Z variables" or "generic variables"), reached through R's Surv()-style formula convention (from the survival package – not a dependency of this one, just the same idea) instead of a separate argument.

Usage

lvl(...)

Arguments

Details

lvl() is only meaningful inside a tda_qreg formula – it is not a real transformation, just a marker tda_qreg looks for and removes before fitting, the same way Surv() is not itself computed by coxph(), only recognised by it.

Value

A matrix, for model.matrix's benefit – not meaningful on its own.

Examples

# Daganzo's classic mode-choice example (examples/exam/qr4.cf in TDA's
# own sources): commute time by car (Z1), bus (Z2), and train (Z3), one
# coefficient for time itself, shared across all three modes
# (seeded synthetic version: choices generated from utility
# -0.3 * time + Gumbel noise, so the estimate recovers about -0.3)
set.seed(21)
n <- 80
d <- data.frame(Z1 = runif(n, 5, 30), Z2 = runif(n, 5, 30),
                Z3 = runif(n, 5, 30))
u <- -0.3 * as.matrix(d) + matrix(-log(-log(runif(3 * n))), n)
d$Y <- max.col(u)
tda_qreg(Y ~ lvl(Z1, Z2, Z3), d, model = "multinomial_logit",
         nq = 3, intercept = FALSE)$estimates

Description

TDA's graph-drawing command, pltree, only draws trees. For anything else – a cycle, more than one component, or a directed graph – this lays the nodes out on a circle and draws it with tda_pl_graph instead, which takes explicit coordinates and so can draw any graph.

Usage

## S3 method for class 'tda_graph'
plot(
  x,
  width = 110,
  height = 80,
  file = "graph.ps",
  newpage = TRUE,
  layout = c("auto", "circle", "tree"),
  ...
)

Arguments

Value

The parsed drawing operations, invisibly.

See Also

Other graph analysis: tda_dmet(), tda_g(), tda_g_analyses, tda_graph(), tda_locate_line(), tda_ptree(), tda_ucl()

Examples

e <- data.frame(from = c(1, 1, 2), to = c(2, 3, 4), value = 1)
g <- tda_graph(e, directed = FALSE)  # 4 nodes, 3 edges: already a tree
pdf(NULL); plot(g); dev.off()

# a graph with a cycle is not a tree, and this one is directed besides --
# both are fine with the default circular layout
e2 <- data.frame(from = c(1, 1, 2, 2, 3, 4), to = c(2, 3, 3, 4, 4, 5))
g2 <- tda_graph(e2, directed = TRUE)
pdf(NULL); plot(g2); dev.off()

Description

Base graphics, drawn from the tables the estimators already return, so no plotting dependency is needed. Censored times are marked, a band is drawn where a standard error is available, and groups get a legend.

Usage

## S3 method for class 'tda_ple'
plot(
  x,
  conf.int = TRUE,
  mark.censored = TRUE,
  col = NULL,
  lty = 1,
  lwd = 1.5,
  xlab = "Time",
  ylab = "Survivor function",
  main = NULL,
  legend = TRUE,
  ylim = c(0, 1),
  ...
)

## S3 method for class 'tda_ltb'
plot(x, ...)

## S3 method for class 'tda_rate'
plot(
  x,
  what = c("survivor", "rate", "density"),
  xlab = "Time",
  ylab = NULL,
  col = "black",
  lwd = 1.5,
  main = NULL,
  ...
)

Arguments

Value

The object, invisibly.

See Also

Other plotting: tda_check_ps(), tda_pl(), tda_pl_arc(), tda_pl_axis(), tda_pl_graph(), tda_pl_hist, tda_pl_panel(), tda_pl_regression(), tda_pl_scatter(), tda_plot_ps(), tda_ps(), tda_ps3(), tda_read_ps()

Examples

set.seed(40)
d <- data.frame(t = round(rexp(60, 0.1), 1) + 0.5, s = rbinom(60, 1, 0.8))
km <- tda_km(Surv(t, s) ~ 1, d)
plot(km)

Description

sdplot needs a PostScript coordinate system as well as the spatial data, so the ranges have to be given.

Usage

## S3 method for class 'tda_spatial'
plot(
  x,
  xlim = NULL,
  ylim = NULL,
  width = 100,
  height = NULL,
  newpage = TRUE,
  what = c("objects", "hull", "hull_per_polygon"),
  symbol = NULL,
  size = NULL,
  lty = NULL,
  lw = NULL,
  gray = NULL,
  ...
)

Arguments

Details

Points need a marker or nothing is drawn: give symbol (or another TDA symbol number). This is the same trap noted for tda_sd_plot3 – without a marker or a zvar attribute, sdplot draws nothing at all rather than erroring.

Value

The parsed drawing operations, invisibly.

See Also

Other spatial analysis: tda_map(), tda_polygons(), tda_read_dbf(), tda_read_shapefile(), tda_read_spatial, tda_sd(), tda_sd_analyses, tda_sd_ps3(), tda_spatial()

Examples

d <- data.frame(id = 1:4, x = c(0, 1, 1, 0), y = c(0, 0, 1, 1))
s <- tda_spatial(d)
plot(s, symbol = 1)  # a marker is needed; points draw nothing without one

Description

TDA's arcd: opens access described by an archive description file, or closes it with file = NULL.

Usage

tda_arcd(file = NULL, ...)

Arguments

Value

the printed output, invisibly.

Examples

# an archive description file names the zoo archive and then one
# line per member: number, name, type (1 data, 2 variable
# descriptions), record length, records, variables
dr <- tempfile("tda"); dir.create(dr)
zoo <- system.file("extdata", "tda.zoo", package = "tdaR")
file.copy(zoo, file.path(dr, "tda.zoo"))
writeLines(c("tda.zoo", "1 adata.dat 1 24 20 3", "2 avar.dat 2 40 3 0"),
           file.path(dr, "arc.ad"))
out <- tda_arcd("arc.ad", dir = dr)
cat(grep("archive|Checking|adata", out, value = TRUE), sep = "\n")

Description

TDA's arcv on an open archive. It requires an open archive (see tda_arcd) containing a type-2 (variable description) member.

Usage

tda_arcv(out, data_file = NULL, archive = NULL, ...)

Arguments

Value

the printed output, invisibly.

Examples

d <- tempfile(); dir.create(d)
for (f in c("tda.zad", "tda.zoo"))
    file.copy(system.file("extdata", f, package = "tdaR"), d)
tda_arcv("vars.vd", archive = "tda.zad", dir = d)
readLines(file.path(d, "vars.vd"))

Description

TDA's arcvc. It requires an open archive (see tda_arcd) containing a type-2 (variable description) member.

Usage

tda_arcvc(file, rewrite = NULL, suffix = c("number", "lower", "upper"), ...)

Arguments

Value

the printed output, invisibly.

Examples

# the file arcv writes, checked for duplicate variable names
d <- tempfile(); dir.create(d)
for (f in c("tda.zad", "tda.zoo"))
    file.copy(system.file("extdata", f, package = "tdaR"), d)
tda_arcv("vars.vd", archive = "tda.zad", dir = d)
# arcv writes a template with every definition commented out, for the
# caller to uncomment; arcvc reports its verdict on it either way
out <- tda_arcvc("vars.vd", dir = d)
cat(grep("variable|once", out, value = TRUE), sep = "\n")

Description

TDA's bfc: Y as a Boolean function of binary X's.

Usage

tda_bfc(y, x, undefined = c("dontcare", "as_one", "as_zero"), ...)

Arguments

Value

the printed output (implicants).

Examples

x <- expand.grid(a = 0:1, b = 0:1)
r <- tda_bfc(y = as.integer(x$a & x$b), x = x)
cat(head(tda_payload(r), 12), sep = "\n")

Description

Which combinations of binary conditions produce the outcome? TDA's bfa builds the truth table of y over the condition variables, finds all prime implicants (Quine-McCluskey), and selects minimal covers. This is the machinery behind crisp-set QCA-style analyses.

Usage

tda_boolean_min(
  y,
  conditions,
  data = NULL,
  undefined = c("dont_care", "true", "false"),
  algorithm = c("lawler", "petrick", "lawler2", "none"),
  max_implicants = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Value

A list: $selections, a data frame with one row per term of each minimal cover – selection number and one column per condition holding 1, 0, or NA for "does not matter" – and $expression, the printed symbolic forms.

See Also

Other optimization: tda_mlp(), tda_mlpi(), tda_mls(), tda_mqap(), tda_mqp()

Examples

d <- expand.grid(X1 = 0:1, X2 = 0:1, X3 = 0:1)
d$Y <- as.integer(d$X1 & d$X2 | d$X3)
tda_boolean_min("Y", c("X1", "X2", "X3"), data = d)$selections

Description

An interval regression identifies its coefficient only up to an interval, so coef() would be misleading. This returns the bounds TDA reports.

Usage

tda_bounds(x, what = c("beta", "alpha", "mean_y", "mean_x", "var_x"))

Arguments

Value

A named vector with lower and upper, or NULL.

See Also

Other interval-valued data: tda_idf(), tda_ilsreg(), tda_imean(), tda_imreg(), tda_inpreg(), tda_ivar1(), tda_ivreg(), tda_sddf()

Examples

set.seed(37)
n <- 40
xlo <- round(rnorm(n), 2); xhi <- xlo + round(runif(n, 0.5, 1.5), 2)
ylo <- round(1 + 0.6 * xlo + rnorm(n, sd = 0.3), 2)
yhi <- ylo + round(runif(n, 0.5, 1.5), 2)
d <- data.frame(xlo = xlo, xhi = xhi, ylo = ylo, yhi = yhi)
ivf <- tda_ivreg(iv(ylo, yhi) ~ iv(xlo, xhi), d)
tda_bounds(ivf, "beta")
tda_bounds(ivf, "mean_x")

# ilsreg: an interval response, a point regressor -- beta and alpha (the
# intercept) both come back bounded, not just the slope
set.seed(2)
d2 <- data.frame(x = round(rnorm(30), 2))
d2$ylo <- 2 + 0.5 * d2$x + rnorm(30, sd = 0.2)
d2$yhi <- d2$ylo + round(runif(30, 0.3, 1), 2)
lsf <- tda_ilsreg(iv(ylo, yhi) ~ x, d2)
tda_bounds(lsf, "beta")
tda_bounds(lsf, "alpha")

Description

brr: builds the replicate-weight matrix balanced repeated replication needs – a variance-estimation technique from survey statistics that works by refitting an estimator on nr half-sample replicates of the data (nr, the number of replications required, is determined from strata and secu and reported back) and taking the spread across them as the variance, rather than a model-based formula. This returns the design itself – TDA's orthogonal (Hadamard-based) construction, not practical to reproduce independently in R – so it can be applied to whatever estimator is being replicated; actually refitting a model nr times and computing the BRR variance from the results is ordinary R code once the design is in hand, not shown here since it depends entirely on what is being estimated.

Usage

tda_brr(strata, secu, dir = tempfile("tda"), ...)

Arguments

Details

Not every (strata, secu) combination succeeds: secu must be at least 2 (one stratum, one unit gives nothing to replicate), and at least one combination that satisfies every other documented constraint (brr = 4,2) can still fail with “no success in calculating orthogonal weights”: 2, 3, and 8 strata at secu = 2 all succeed, 4 does not. This is a limitation in TDA's Galois-field construction for specific sizes, not something this wrapper can lift.

Value

A list: replications (nr, how many are needed), and design, the strata by replications matrix itself (TDA's description: design[i, j] is which sampling unit of stratum i to use in replication j).

See Also

Other descriptive statistics: tda_cov(), tda_dstat(), tda_freq(), tda_independence(), tda_ineq(), tda_loglin(), tda_quant(), tda_rcorr(), tda_segr(), tda_subm()

Examples

b <- tda_brr(3, 2)
b$replications
b$design

Description

Small file-diagnostic utilities from TDA's toolkit, from an era before assuming wc, od or split were available: tda_ccnt a character frequency table, tda_lcnt a record-length frequency table, tda_dump a hex dump, and tda_dsplit splits a file into fixed-size parts on disk.

Usage

tda_ccnt(file, dir = tempfile("tda"), ...)

tda_lcnt(file, noc = NULL, dir = tempfile("tda"), ...)

tda_dump(file, nc = NULL, offset = NULL, dir = tempfile("tda"), ...)

tda_dsplit(file, len = NULL, dir = tempfile("tda"), ...)

Arguments

Value

tda_ccnt and tda_lcnt return a data frame. tda_dump returns the hex dump as a character vector, one line each. tda_dsplit returns the paths of the parts it wrote.

See Also

Other TDA infrastructure: tda_help(), tda_output(), tda_read_table(), tda_run(), tda_write_data()

Examples

f <- tempfile()
writeLines(c("hello world", "foo bar baz", "x"), f)
tda_ccnt(f)
tda_lcnt(f)
cat(tda_dump(f), sep = "\n")

f2 <- tempfile()
writeLines(rep("0123456789", 10), f2)
tda_dsplit(f2, len = 20)

Description

Counts the drawing operators in the file and the operations the parser produced, and reports anything unaccounted for. Reading a plot back is easy to get half right – an operator that is silently ignored looks exactly like one that is not there – so this compares the two directly rather than leaving it to the eye.

Usage

tda_check_ps(file, ...)

Arguments

Value

A data frame with one row per operator: how often it appears in the file, whether the parser handles it, and how many operations resulted.

See Also

Other plotting: plot.tda_ple(), tda_pl(), tda_pl_arc(), tda_pl_axis(), tda_pl_graph(), tda_pl_hist, tda_pl_panel(), tda_pl_regression(), tda_pl_scatter(), tda_plot_ps(), tda_ps(), tda_ps3(), tda_read_ps()

Examples

p <- tda_ps(xlim = c(0, 10), ylim = c(0, 10))
p <- tda_pl_circle(p, at = c(5, 5), r = 2)
tda_check_ps(tda_ps_file(p))

Description

TDA's clpyr on a distance matrix.

Usage

tda_clpyr(d, linkage = c("single", "complete"), ...)

Arguments

Value

the printed output.

Examples

r <- tda_clpyr(as.matrix(dist(c(0, 1, 5, 6))))
# the top of the pyramid joins everything: 4 objects spanning 1..4
tail(r$pyramid, 1)

Description

TDA's clustering commands, all of which take a dissimilarity and so accept a dist – including the one tda_seqm returns.

Usage

TDA_CLUSTER

tda_cluster(
  d,
  method = "hierarchical",
  algorithm = c("centers", "diameter"),
  min_size = NULL,
  max_splits = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

The default, "hierarchical" (hcld), writes a ragged table – a node, its member count, then that many member indices – which is not a rectangular data.frame, so it is parsed separately into clusters, a named list of member-index vectors, one per node of the divisive tree (the root first, then each split). "unidimensional" and "additive" report their result only as text (a final partition or a fitted tree), never through df=; read $run$output for those.

Value

An object carrying a table, or for "hierarchical" clusters instead (see Details). For "hierarchical" the object also carries merges: one row per split, the two clusters, their sizes and their diameters.

See Also

Other clustering: TDA_MDS, tda_conjoint(), tda_cutree(), tda_dma(), tda_pdatd()

Examples

d <- dist(matrix(c(1, 2, 3, 8, 9, 10), ncol = 1))
tda_cluster(d)$clusters  # the divisive tree: root, then each split

Description

TDA's com: enumerate tuples, sets, permutations or partitions into an output file.

Usage

tda_com(
  pattern = c("tuples", "subsets", "m_tuples", "permutations", "partitions",
    "partitions_m", "bipartitions"),
  n = 1,
  m = 1,
  ...
)

Arguments

Value

the output file's lines.

Examples

# the six 2-element subsets of {0, 1, 2, 3}: a row number, then the members
tda_com("subsets", n = 4, m = 2)
# the five partitions of 4: 4, 3+1, 2+2, 2+1+1, 1+1+1+1
tda_com("partitions", n = 4)

Description

tda_combine_ps reaches TDA's dplot command, which lays out a grid of previously-created PostScript files (each from its own tda_ps() session, already run) into a single output file – one row per element of rows, side by side within a row. Unlike every other tda_pl_* function here, this is not a command added to an in-progress session: dplot combines files that already exist on disk into a new one of its own, so this takes file paths (or tda_ps objects that have already been run, via tda_ps_file) and runs standalone, the way tda_run does.

Usage

tda_combine_ps(
  rows,
  width = 120,
  height = 80,
  origin = c(10, 10),
  file = "combined.ps",
  ...
)

Arguments

Value

The result of tda_run, whose own run directory holds file.

Examples

# two small plots, then both side by side in one file
p1 <- tda_ps(xlim = c(0, 10), ylim = c(-1, 1), width = 55, height = 40)
p1 <- tda_pl_axes(tda_pl_function(p1, "sin(x)"))
p2 <- tda_ps(xlim = c(0, 10), ylim = c(0, 100), width = 55, height = 40)
p2 <- tda_pl_axes(tda_pl_function(p2, "x * x"))
r <- tda_combine_ps(list(c(tda_ps_file(p1), tda_ps_file(p2))))
file.exists(file.path(r$dir, "combined.ps"))

Description

TDA's conj (alternating least squares) on a preference variable and factor columns.

Usage

tda_conj(y, x, algorithm = c("monotone", "lp"), ...)

Arguments

Value

the printed output.

Examples

x <- expand.grid(a = 1:2, b = 1:2)
r <- tda_conj(y = c(4, 3, 2, 1), x = x)
cat(tail(tda_payload(r), 8), sep = "\n")   # the part-worths

Description

nmca is conjoint measurement: it takes a rank order over profiles and the categorical factors describing them, and estimates the part-worth of each factor level. It draws a Hasse diagram of the order as part of the output.

Usage

tda_conjoint(formula, data, options = list(), dir = tempfile("tda"), ...)

Arguments

Value

An object carrying the run; the diagram and the estimates are in $run$output.

See Also

Other clustering: TDA_CLUSTER, TDA_MDS, tda_cutree(), tda_dma(), tda_pdatd()

Examples

d <- expand.grid(price = 1:2, brand = 1:2)
d$rank <- c(4, 2, 3, 1)
cj <- tda_conjoint(rank ~ price + brand, d)
cat(tda_payload(cj$run), sep = "\n")   # the nmca section itself

Description

tda_rate takes constraints in TDA's notation, where parameters are numbered in the order they are estimated: "b3 - b10 = 0" makes the third and tenth equal. This builds such a string from a fitted model's coefficient names, so the numbering does not have to be counted out by hand.

Usage

tda_constrain(fit, ...)

Arguments

Value

A character vector of constraints for tda_rate.

See Also

Other rate models: tda_control(), tda_dple(), tda_rates(), tda_split(), tda_survivor(), tda_transitions(), vcov.tda_fit()

Examples

d <- tda_rrdat(states = 4)
f <- tda_rate(Surv(TFP, DES) ~ EDU + PRES, d)
tda_constrain(f, c("0->1: EDU", "0->2: EDU"))

Description

Convergence control

Usage

tda_control(
  algorithm = NULL,
  tolsg = NULL,
  tolsp = NULL,
  tolg = NULL,
  tolf = NULL,
  maxit = NULL,
  tol_param = NULL,
  tol_var = NULL,
  tol_reduction = NULL,
  step_length = NULL,
  step_reduction = NULL,
  step_min = NULL,
  maxit_line = NULL,
  criterion = NULL,
  derivatives = NULL
)

tda_strict(tol = 1e-12)

tda_converged(x)

Arguments

Value

A list of options for control=; tda_converged returns TRUE, FALSE, or NA where TDA reports nothing.

Criteria

TDA's default minimiser stops on the norm of the gradient against an absolute tolerance, plus a loose disjunction that also stops when the change in function value is small. Algorithms 7 and 8 switch to optstp(), the Dennis-Schnabel test on the largest relative gradient and relative parameter step. That is scale-invariant, closer to what R's optimisers do, and gives results that reproduce across platforms. It is reachable only through algorithm: TDA's crit= accepts 1 to 3 and cannot select it.

The default is TDA's, so a plain call reproduces a plain command file. tda_strict() opts in.

A caution

Algorithms 7 and 8 also lower the iteration cap to 20 when derivatives are numerical, against 100 otherwise, so asking for a tighter tolerance can turn a converged fit into a truncated one. Every fit reports how it stopped, warns if it did not converge, and tda_converged() reads that back.

See Also

Other rate models: tda_constrain(), tda_dple(), tda_rates(), tda_split(), tda_survivor(), tda_transitions(), vcov.tda_fit()

Examples

set.seed(36)
d <- data.frame(x = rnorm(200))
d$y <- rbinom(200, 1, plogis(0.3 + 0.8 * d$x))
f <- tda_qreg(y ~ x, d, control = tda_control(maxit = 50, tolg = 1e-8))
coef(f)
tda_converged(f)

# tda_strict: the same idea, one shared tolerance for both algorithm 7/8
# criteria at once
f2 <- tda_qreg(y ~ x, d, control = tda_strict())
coef(f2)

Description

The $n - 1$-denominator covariance matrix, as stats::cov() computes it, extended to case weights when a cwt command is active. tda_corr is its correlation counterpart.

Usage

tda_cov(..., options = list(), dir = tempfile("tda"))

Arguments

Value

An object carrying matrix, the covariance matrix with the variables' names on both dimensions.

See Also

Other descriptive statistics: tda_brr(), tda_dstat(), tda_freq(), tda_independence(), tda_ineq(), tda_loglin(), tda_quant(), tda_rcorr(), tda_segr(), tda_subm()

Examples

set.seed(26)
d <- data.frame(x = rnorm(40), y = rnorm(40))
d$y <- d$y + 0.6 * d$x
tda_cov(d)$matrix
cov(d)  # matches

Description

Builds a dendrogram with TDA's hcls and cuts it with hclsp.

Usage

tda_cutree(
  d,
  nlev = 2,
  nodes = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

Unlike R's cutree, which returns one group label per observation, hclsp reports the cut as a set of nodes: table has one row per resulting cluster, its TDA node number, how many leaves (original observations) it holds, and its parent node – not a per-observation vector. Getting a label for each observation means walking merge (the dendrogram hcls built, one row per step) down from a node in table to the leaves under it; this wrapper does not do that walk for you.

Value

An object carrying merge, the dendrogram, and table, one row per cluster in the cut (node, leaves, parent).

See Also

Other clustering: TDA_CLUSTER, TDA_MDS, tda_conjoint(), tda_dma(), tda_pdatd()

Examples

set.seed(35)
pts <- rbind(matrix(rnorm(10, 0), 5, 2), matrix(rnorm(10, 8), 5, 2))
ct <- tda_cutree(dist(pts), nlev = 2)
ct$table    # two clusters of five leaves each
ct$merge    # the nine merge steps that built the dendrogram

Description

tda_cwt_norm reproduces TDA's cwt(wnorm=s)=W; (or the bare-flag cwt(wnorm)=W;) as a plain vector transform, for weights= on tda_qreg, tda_glm, tda_ple, tda_ltb and tda_rate – none of which expose wnorm as an argument, since this is all it does: fitted coefficients come back identical whether the weights are rescaled this way first or not, only the standard errors change (rescaling is a constant multiplier on the whole weighted log-likelihood, so it moves the curvature at the optimum, not where the optimum is). cwt(wnorm=s)=W; on the raw weights and plain cwt=W; on tda_cwt_norm(W, s) give the same fit, to the last printed digit.

Usage

tda_cwt_norm(w, s = length(w))

Arguments

Value

w, rescaled so it sums to s.

Examples

w <- c(2, 2, 2, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10)  # sums to 74
tda_cwt_norm(w)         # cwt(wnorm)=W;    -- rescaled to sum 13 (length(w))
tda_cwt_norm(w, 100)    # cwt(wnorm=100)=W;

Description

TDA's standalone dblock on the current-run data: blocks are maximal runs of identical values. Mostly useful inside longer tda_run pipelines; for block-mode variable creation put dblock= inside nvar (see the package tests).

Usage

tda_dblock(data, by, commands = character(), ...)

Arguments

Value

the printed output, invisibly.

Examples

out <- tda_dblock(data.frame(G = c(1, 1, 2), X = 1:3), by = "G",
                  commands = "dstat = X;")
# blocks are defined by G; dstat runs under the block structure
grep("^X", out, value = TRUE)

Description

Runs nvar on a data frame with new variables defined by TDA expressions, and returns the frame with them added. This is the way to the operators of the manual's 5.2.6 that have no R counterpart: the change and count operators (change, cntch), the aggregates over blocks (gcnt, grec, gmean, gstd, ...), TDA's random numbers, and so on.

Usage

tda_derive(data, ..., block = NULL, options = list(), dir = tempfile("tda"))

Arguments

Value

data with the new columns, from the export channel.

Examples

d <- data.frame(ID = c(1, 1, 2), S = c(1, 3, 2))
tda_derive(d, CH = "change(S)", CN = "cntch(S)")
tda_derive(d, CN = "cntch(S)", block = "ID")

Description

Every error, warning, note or "cannot" line TDA printed, on stdout or stderr, for a run or a fitted model. TDA reports problems by printing them and carrying on, so a fit can come back with a full table and a message about what it could not do; nothing here is dropped from tda_output(), this is the short list.

Usage

tda_diagnostics(x)

Arguments

Value

a character vector, empty when TDA reported nothing.

Examples

# a clean run reports nothing
tda_diagnostics(tda_run("nvar(noc = 5, X = case);"))

# dstat asks for a variable that was never defined: TDA prints the
# error and continues, so the run returns and the message is here
r <- tda_run(c("nvar(noc = 5, X = case);", "dstat = Y;"))
tda_diagnostics(r)

# a fit stopped after one iteration still returns a coefficient
# table; the problem is reported alongside it
d <- data.frame(x = c(-1, -0.5, 0, 0.5, 1, 1.5), y = c(0, 0, 1, 0, 1, 1))
f <- suppressWarnings(tda_qreg(y ~ x, d, control = tda_control(maxit = 1)))
tda_diagnostics(f)

Description

dma: seven related methods for reducing a data matrix to a handful of derived scores, selected with alg – principal components (1, from the raw data; 2, from a covariance or correlation matrix already computed), factor analysis (3, 4, the same two starting points but with eigenvectors rescaled as loadings), dual scaling and correspondence analysis (5, 6, for a frequency table rather than a data matrix), and a direct SVD-based projection (7). Verified against prcomp for algorithm 1: eigenvalues match exactly, eigenvectors match up to sign (an ordinary, harmless ambiguity in any eigendecomposition).

Usage

tda_dma(
  data,
  alg = 1,
  preprocess = NULL,
  ns = NULL,
  variables = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

TDA's v= (choosing a subset of variables from a larger data matrix) needs variable names that start uppercase, the same rule as anywhere else in TDA – lowercase R column names are displayed correctly when declared but cannot be referenced again later in the same script, including by v= itself, failing with a “Syntax error or undefined variables” that gives no hint the name's case is the problem. variables handles this automatically; there is no need to rename columns yourself first.

Value

A list with values (the eigenvalues), percent (each eigenvalue's share of the total, as TDA prints it), vectors (the eigenvectors or factor loadings, one column per component, row names taken from data), and scores (the derived score for each case, one column per component) – scores is empty for algorithms where TDA does not calculate it (2 and 4).

See Also

Other clustering: TDA_CLUSTER, TDA_MDS, tda_conjoint(), tda_cutree(), tda_pdatd()

Examples

set.seed(1)
d <- data.frame(x1 = rnorm(30), x2 = rnorm(30))
d$x3 <- d$x1 + d$x2 + rnorm(30, sd = 0.3)   # correlated with both
pc <- tda_dma(d)
pc$values                 # most of the variance is in the first PC
pc$vectors                # loadings: x3 dominates it
head(pc$scores)

# a subset of the variables, instead of all of them
tda_dma(d, variables = c("x1", "x2"))$values

Description

dmet: an undirected valued graph's edge values, read as a distance matrix, do not always satisfy the triangle inequality (a direct edge can be "longer" than going via a third node). This finds the smallest constant that, added to every distance, makes the whole matrix satisfy it – $c = \max(0, \max_{i,j,k} d_{ij} - d_{ik} - d_{kj})$ – and returns both the original and the corrected matrix.

Usage

tda_dmet(g, options = list(), dir = tempfile("tda"), ...)

Arguments

Value

A list: original and modified, both distance matrices, and constant, the value added to get from one to the other.

See Also

Other graph analysis: plot.tda_graph(), tda_g(), tda_g_analyses, tda_graph(), tda_locate_line(), tda_ptree(), tda_ucl()

Examples

e <- data.frame(from = c(1, 2, 3), to = c(2, 3, 1), value = c(1, 2, 5))
g <- tda_graph(e, directed = FALSE)   # 1-3 direct (5) exceeds
                                      # 1-2-3 via node 2 (1+2=3)
tda_dmet(g)

Description

TDA's dmet1: adjusts a distance matrix (given as an undirected valued graph) toward metricity.

Usage

tda_dmet1(d, tolerance = 1e-04, max_iterations = 20, ...)

Arguments

Value

the output file's lines.

Examples

# 5 breaks the triangle inequality (1 + 2 < 5), so the
# metricized matrix shortens it to 3
d <- rbind(c(0, 1, 5), c(1, 0, 2), c(5, 2, 0))
r <- tda_dmet1(d)
r$modified

Description

tda_ple is the product-limit estimator, which is Kaplan-Meier; tda_km is an alias for it. tda_ltb groups into intervals. tda_dple and tda_dltb are the discrete time versions, and tda_diple handles interval censored data.

Usage

tda_dple(formula, data, dir = tempfile("tda"), ...)

tda_dltb(formula, data, dir = tempfile("tda"), ...)

tda_diple(data, start, lower, upper, status, dir = tempfile("tda"), ...)

tda_ltb(
  formula,
  data,
  tp = NULL,
  cfrac = NULL,
  weights = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

tda_ple(
  formula,
  data,
  compare = FALSE,
  quantiles = NULL,
  at_survival = NULL,
  weights = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

tda_km(formula, data, ...)

Arguments

Details

The survivor function's standard errors are the usual Greenwood estimate (the manual gives the formula in section 6.5.2), the same one survival::survfit reports, so the two agree on identical data.

Value

An object carrying table, blocks with one entry per group or transition, and for the product limit median. Be precise about the pair: table is deliberately only the first block, the convenience view for the common single-group call; a grouped or multi-transition fit's full result is blocks (named by group), and anything that must see everything – comparisons, exports, your summaries – should read do.call(rbind, fit$blocks), never table alone. The product-limit columns are time, events, censored, n.risk, survivor, std.err and cum.rate.

Censoring and competing risks

The second argument of Surv() is the destination state, not simply an event indicator. Zero means no transition, so it is a censored case. A value above 1 is a second destination, which makes the data competing risks, and the estimate is returned per transition in $blocks.

See Also

Other rate models: tda_constrain(), tda_control(), tda_rates(), tda_split(), tda_survivor(), tda_transitions(), vcov.tda_fit()

Examples

d <- data.frame(t = c(4, 3, 1, 1, 2, 2, 3, 5), s = c(1, 1, 1, 0, 1, 1, 0, 1))
km <- tda_km(Surv(t, s) ~ 1, d)
km$table
plot(km)

# ltb: the same data grouped into intervals rather than at each event time
tda_ltb(Surv(t, s) ~ 1, d, tp = seq(0, 6, 2))$table

# ple: tda_km's underlying function, with two features tda_km does
# not expose -- compare = TRUE tests whether group survivor functions
# differ, and quantiles reads the curve at chosen times
set.seed(1)
d1b <- data.frame(t = round(rexp(40, 0.1), 1) + 0.5, s = rbinom(40, 1, 0.8),
                  grp = rep(c("A", "B"), 20))
pf <- tda_ple(Surv(t, s) ~ factor(grp), d1b, compare = TRUE,
             quantiles = c(5, 10))
pf$comparison  # a log-rank test between the two groups' curves
pf$quantiles   # the survivor function at t = 5 and t = 10, per group

# dple, dltb: the discrete-time versions, for integer time already grouped
# into periods -- t itself is the period here, not a continuous duration
set.seed(44)
d2 <- data.frame(t = sample(1:10, 60, TRUE), s = rbinom(60, 1, 0.8))
tda_dple(Surv(t, s) ~ 1, d2)$table
tda_dltb(Surv(t, s) ~ 1, d2)$table

# diple: the event is known only to fall between lower and upper
set.seed(45)
d3 <- data.frame(start = 0, lower = sample(1:8, 40, TRUE))
d3$upper <- d3$lower + sample(1:3, 40, TRUE)
d3$status <- rbinom(40, 1, 0.85)
tda_diple(d3, start = "start", lower = "lower", upper = "upper",
         status = "status")$table

Description

These take data rather than a formula, following cor and summary: there is no response and no model here. A data frame, a matrix, or loose vectors are all accepted, and non-numeric columns are dropped.

Usage

tda_dstat(..., by = NULL, options = list(), dir = tempfile("tda"))

tda_corr(..., options = list(), dir = tempfile("tda"))

Arguments

Details

tda_dstat reports each variable's minimum, maximum, mean, standard deviation (denominator $n - 1$, or the weight total minus one under case weights) and sum. tda_corr and tda_cov are the ordinary product-moment correlation and the $n - 1$-denominator covariance – cor() and cov() – extended to case weights when a cwt command is active (manual, section 6.2.6). All of them compute over every case given, including any -5 standing in for an NA – see ?tdaR on missing values.

Value

An object carrying a table, or for tda_corr a matrix.

See Also

Other descriptive statistics: tda_brr(), tda_cov(), tda_freq(), tda_independence(), tda_ineq(), tda_loglin(), tda_quant(), tda_rcorr(), tda_segr(), tda_subm()

Examples

set.seed(1)
d <- data.frame(x = round(rnorm(40, 10, 3), 1))
d$y <- round(2 + 0.8 * d$x + rnorm(40, sd = 2), 1)
tda_dstat(d)
tda_corr(d$x, d$y)

Description

Runs TDA's etest, which computes the eigenvalues and eigenvectors of a square matrix and reports the residual of each. It is a test of TDA's eigen routines rather than a general facility – eigen is what you want for ordinary work – but it is the way to check what TDA itself computes for a matrix a model is about to be fitted on.

Usage

tda_eigen(
  x,
  algorithm = 1,
  max_iter = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

The eigenvalues come from an exporter rather than from the printed table, so they arrive at full precision instead of rounded to the print format.

Value

A data frame with one row per eigenvalue: re and im, its real and imaginary parts. attr(x, "run") carries the run, whose output holds the eigenvectors and residuals.

Examples

m <- diag(c(2, 3, 5))
tda_eigen(m)

Description

tda_episodes returns the episode data edef builds from a formula, which is what every model in this package is actually fitted to. It is the way to check that a Surv() specification means what you intended before trusting a fit.

Usage

tda_episodes(
  formula,
  data,
  vars = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

tda_state_dist(
  formula,
  data,
  times,
  id = NULL,
  spell = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

tda_state_dist gives the distribution over states at a set of time points – how many cases occupy each state when – which is the natural summary of multi-state data.

Value

A data frame of episodes, or of state proportions by time.

See Also

Other episodes: tda_esort()

Examples

d <- data.frame(t = c(4, 3, 1, 5), s = c(1, 0, 1, 1), x = c(2, 1, 3, 1))
tda_episodes(Surv(t, s) ~ x, d)
tda_state_dist(Surv(t, s) ~ x, d, times = c(1, 3, 5))

Description

TDA was written when a data file could be larger than the machine's memory, and these commands work a record at a time on an ASCII file rather than loading it: esort sorts the records, eskip drops character positions, eselect keeps the records whose key appears in a second file, emerge joins sorted files on a key. They are file-in, file-out, on raw text records addressed by character position, and they are that here too: name the input and the output file and the output path comes back, invisibly, with the run as its "run" attribute. Nothing is read back; for a data frame, use tda_read_table on the file you asked for.

Usage

tda_esort(file, keys, out, options = list(), dir = tempfile("tda"))

tda_eskip(file, drop, out, options = list(), dir = tempfile("tda"))

tda_eselect(
  file,
  with,
  keys,
  with_keys = keys,
  out,
  options = list(),
  dir = tempfile("tda")
)

tda_emerge(
  file,
  with,
  keys,
  with_keys = keys,
  out,
  options = list(),
  dir = tempfile("tda")
)

tda_ejoin(
  data,
  id,
  start,
  end,
  state,
  with = NULL,
  with_id = NULL,
  with_start = NULL,
  with_end = NULL,
  with_state = NULL,
  options = list(),
  dir = tempfile("tda")
)

Arguments

Details

tda_emerge is a merge join and requires both inputs sorted in ascending order on the key; it says so and stops otherwise, naming the two records that are out of order.

tda_ejoin is different in kind: it joins episode data into levels, takes data frames, and returns one through the export channel.

Value

The output path, invisibly, for the four file utilities; a data frame for tda_ejoin.

See Also

Other episodes: tda_episodes()

Examples

f <- tempfile(); writeLines(c("103 30", "101 10", "102 20"), f)
o <- tempfile()
tda_esort(f, keys = c(1, 3), out = o)
readLines(o)
tda_eskip(f, drop = c(1, 4), out = o)
readLines(o)

# eselect/emerge need their inputs already sorted on the key
f1 <- tempfile(); writeLines(c("101 10", "102 20", "103 30"), f1)
f2 <- tempfile(); writeLines(c("101 1", "103 3"), f2)
tda_eselect(f1, with = f2, keys = c(1, 3), out = o); readLines(o)
tda_emerge(f1, with = f2, keys = c(1, 3), out = o); readLines(o)

# ejoin: joins one or two episode datasets, splitting overlapping
# spells into levels -- an activity spanning 0-10, and a state that
# changes mid-activity at t=4, joined so each row is a period where
# neither changes
act <- data.frame(id = 1, start = 0, end = 10, state = 1)
st <- data.frame(id = 1, start = c(0, 4), end = c(4, 10), state = c(1, 2))
tda_ejoin(act, id = "id", start = "start", end = "end", state = "state",
         with = st)

Description

evalfi: given an expression and an interval domain for each of its arguments, returns an interval guaranteed to contain the expression's true range over that domain – interval arithmetic, not a numeric evaluation at a point. The bound is conservative rather than tight where a variable appears more than once (the standard "dependency problem" of interval arithmetic: each occurrence is widened independently, as if they were unrelated variables that happened to share a domain) – evalfi(x = c(0, 3), "x*x-2*x+1") returns [-5, 10], not the true range [0, 4] that (x-1)^2 takes over [0, 3]. evalfi1 is the same but for the first derivative.

Usage

tda_evalfi(
  expr,
  ...,
  derivative = FALSE,
  options = list(),
  dir = tempfile("tda")
)

Arguments

Details

^ is not a valid operator in an interval expression – write a power as repeated multiplication (x*x, not x^2).

Value

A named numeric vector, c(lower, upper).

See Also

Other sequence analysis: tda_fml(), tda_frml(), tda_seq_info(), tda_seqgc(), tda_seqm(), tda_seqmd(), tda_seqpe(), tda_seqpm()

Examples

# true range of (x-1)^2 over [0,3] is [0,4] -- this is wider, since x
# appears twice and interval arithmetic cannot see they are the same x
tda_evalfi("x*x-2*x+1", x = c(0, 3))
tda_evalfi("x*x-2*x+1", x = c(0, 3), derivative = TRUE)  # 2x-2, exact
                                                         # here: [-2,4]

Description

tda_expand runs TDA's expm, which expands an expression into its fully written-out form – a model formula with interaction shorthand becomes the explicit list of terms. tda_mparse runs mparse, which parses a matrix expression and reports how it read it, without evaluating anything.

Usage

tda_expand(expr, setup = NULL, options = list(), dir = tempfile("tda"), ...)

tda_mparse(
  expr,
  setup = NULL,
  mfmt = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

Both are diagnostics: they answer "what did TDA make of what I wrote", which is otherwise only visible in an error message.

Value

A character vector: the lines TDA printed in response. attr(x, "run") carries the run.

Examples

tda_expand("A*B")
tda_mparse("A+B", setup = c("mdefi(2,2,A);", "mdefi(2,2,B);"))

Description

Fits a generalized linear model, taking R's family objects.

Usage

TDA_FAMILIES

tda_glm(
  formula,
  data,
  family = stats::gaussian,
  link = NULL,
  custom_link = NULL,
  domain = NULL,
  trials = NULL,
  weights = NULL,
  intercept = NULL,
  start = NULL,
  predictions = FALSE,
  protocol = FALSE,
  equality = NULL,
  inequality = NULL,
  control = NULL,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Value

An object of class tda_fit, its fit element (shown by summary()) carrying the deviance and residual degrees of freedom, the rank of the design matrix, and – for the families where TDA prints them – the Pearson statistic and the ML- and deviance-based dispersion estimates. Cases with missing values are dropped before the fit, with a message (see ?tdaR). logLik() is NA: TDA's glm() never prints a log-likelihood itself, only Deviance – and it stays that way here rather than being reconstructed from Deviance, since that reconstruction would be a derived, unverifiable-by-the-user computation on top of what TDA itself actually computes.

Families and links

TDA_FAMILIES maps R's names onto TDA's d= (distribution) and link= numbering, from TDA's manual (which tda_help("glm") cannot show – glm has no entry in TDA's help database, despite being fully implemented). The first link listed for each family is the default, and is the canonical one:

TDA's link numbering runs 1 to 8 in total (identity, log, logit, reciprocal/inverse, probit, cloglog, sqrt, and the inverse-Gaussian's quadratic inverse, 1/mu^2) – every one of the 8 appears above, on whichever family it is offered for through R's family() system. A family TDA does not implement, or a link it does not offer for that family, is rejected with the list of what is available rather than silently replaced by a default.

glm's ab= option – the domain for a user-defined link function, one not among the 8 above, given as its TDA expression on the right-hand side of the command rather than a number – is reachable through custom_link/domain below.

See Also

Other regression: TDA_QRMODELS, tda_freg(), tda_gdf(), tda_l1reg(), tda_lsreg(), tda_mlrc_design(), tda_mreg(), tda_nlreg(), tda_npreg(), tda_zreg()

Examples

set.seed(5)
d <- data.frame(x = round(rnorm(40), 2))
d$y <- rbinom(40, 1, plogis(-0.5 + 1.2 * d$x))
fit <- tda_glm(y ~ x, d, family = binomial)
fit
predict(fit)                              # fitted values, on the data
                                           # the model was fit on
predict(fit, newdata = data.frame(x = 1))

# a binomial model of trial data (successes out of a known number of
# trials), not just 0/1 outcomes -- the response is the success count,
# trials the number of trials, not a proportion
set.seed(1)
n <- 60
d2 <- data.frame(x = rnorm(n))
p <- plogis(0.3 + 0.8 * d2$x)
d2$trials <- sample(5:20, n, replace = TRUE)
d2$successes <- rbinom(n, d2$trials, p)
coef(tda_glm(successes ~ x, d2, family = binomial, trials = d2$trials))

# case weights instead (or as well) -- a different thing from
# the trial counts above
d2$imp <- runif(n, 0.5, 2)
coef(tda_glm(successes ~ x, d2, family = binomial, trials = d2$trials,
            weights = "imp"))

# a constrained fit: the two coefficients forced to sum to exactly 1
d3 <- data.frame(x1 = rnorm(60), x2 = rnorm(60))
d3$y <- rbinom(60, 1, plogis(0.6 * d3$x1 + 0.4 * d3$x2))
coef(tda_glm(y ~ x1 + x2, d3, family = binomial,
            equality = "x1 + x2 = 1"))

# family and link as TDA's raw numbers -- useful when translating
# a .cf file directly, which only ever writes them this way. Reaches
# the identical call as family = Gamma(link = "log").
set.seed(2)
d4 <- data.frame(x = rnorm(60))
d4$y <- rgamma(60, shape = 2, rate = 2 / exp(0.5 + 0.3 * d4$x))
coef(tda_glm(y ~ x, d4, family = 4, link = 2))

# a user-defined link -- the logit link, spelled out by hand instead
# of using link = "logit"; reaches the same fit exactly
coef(tda_glm(y ~ x, d, family = binomial,
            custom_link = "log(mue / (1 - mue))", start = c(0, 0)))

Description

fml maximises a log-likelihood you write yourself, which is what TDA offers that a fixed catalogue of models cannot. The likelihood is written in an expression language, evaluated inside TDA rather than as an R function – but as of this version, that expression can be written in ordinary R syntax and translated automatically, rather than TDA's.

Usage

tda_fml(
  definitions,
  data,
  start = NULL,
  control = NULL,
  constraints = NULL,
  residuals = FALSE,
  residual_vars = NULL,
  protocol = FALSE,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Value

An object of class tda_fit. With residuals = TRUE, also carries $residuals (see residuals above).

Writing the likelihood in R syntax

Give definitions as an unevaluated { } block of plain R assignments instead of a character vector, and it is translated to TDA's own expression syntax for you – never actually run as R code (DES, rate and the rest only ever exist inside TDA, so it cannot be). ifelse(cond, yes, no) becomes TDA's if() (TDA's if takes a 2- or 3-argument if(cond, then)/if(cond, then, else) shape matching ifelse exactly, not R's control-flow if), and R's comparison operators (<, >, <=, >=, ==, !=) become TDA's comparison functions (lt, gt, le, ge, eq, ne). Arithmetic (+ - * / ^), &/|, and ordinary function calls (exp, log, negbin, ...) already mean the same thing in both and pass through unchanged.

These six TDA functions are not the strict comparisons their R equivalents are. Every one of them – lt, le, gt, ge, eq, ne alike – carries a fixed tolerance, EPSI2 = 1000 * DBL_EPSILON (about 2.22e-13): le, say, computes x <= y + EPSI2, not x <= y. This can matter on real data: a case in examples/exam/rrdat.1 has PRESN/PRES - 1 equal to 0.2 exactly in exact arithmetic, and the two land on different sides of 0.2 once represented as IEEE 754 doubles – ge(x, 0.2) inside TDA still calls it true, R's own x >= 0.2 on the identical computation does not. Applied consistently across all six, at least, but nowhere in TDA's documentation: not in the manual, not in tda_help()'s operator table, which describes le as plainly as “true if less than, or equal” with no footnote at all. A hand-written R replication of a TDA calculation involving a comparison near a round number (==, or a boundary like 0.2 here) can therefore disagree with TDA's output at that boundary, for a reason with nothing to do with which language computed it – this is a property of TDA's expression evaluator, present whether the comparison was written by hand in TDA's syntax or reached it through this translation.

tda_fml({
    rate = exp(a0 + COHO2 * a1 + COHO3 * a2 + W * a3)
    l1 = ifelse(DES, log(rate), 0)
    fn = l1 - rate * DUR
}, data = d, start = c(-4, 0, 0, 0))

The older, direct form – a character vector of TDA's assignment syntax – still works exactly as before, for anything the translator does not (yet) cover, or for pasting a .cf file's text in directly.

tda_fml(c("xb = b0 + x * b1",
          "fn = -0.5 * (y - xb)^2"), data = d)

estimates b0 and b1 by least squares written as a likelihood. TDA's functions are available, so fn = negbin(alpha, gamma, y) and the rest work.

See Also

Other sequence analysis: tda_evalfi(), tda_frml(), tda_seq_info(), tda_seqgc(), tda_seqm(), tda_seqmd(), tda_seqpe(), tda_seqpm()

Examples

set.seed(1)
d <- data.frame(x = 1:20, y = 2 + 0.5 * (1:20) + rnorm(20))
tda_fml(c("xb = b0 + x * b1", "fn = -0.5 * (y - xb)^2"), d)

# the same model, written in R syntax instead
tda_fml({
    xb = b0 + x * b1
    fn = -0.5 * (y - xb)^2
}, d)

# a compound Poisson (negative binomial) count regression -- TDA has no
# separate command for this (unlike its plain Poisson, which is a glm()
# family), but negbin(alpha, gamma, k) is a built-in log-likelihood
# function in its expression language, reachable through fml like any
# other: gamma is the mean, alpha = 1/sigma the dispersion
set.seed(7)
d2 <- data.frame(x = rnorm(200))
d2$ndi <- rnbinom(200, mu = exp(0.5 + 0.8 * d2$x), size = 2)
coef(tda_fml({
    xb = b0 + x * b1
    gamma = exp(xb)
    sigma = exp(sig)
    alpha = 1 / sigma
    fn = negbin(alpha, gamma, ndi)
}, d2, start = list(b0 = 0, b1 = 0, sig = 0)))

Description

freg fits a regression function you write yourself, in the same way tda_fml takes a likelihood. The expression is TDA's language, not R's.

Usage

tda_freg(
  definitions,
  data,
  start = NULL,
  control = NULL,
  constraints = NULL,
  residuals = FALSE,
  residual_vars = NULL,
  protocol = FALSE,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

tda_nlreg's expr argument also fits a user-written nonlinear function now, using nlreg's dedicated algorithm (proper standard errors, and orthogonal distance regression as an option) rather than a general-purpose minimiser – prefer it for an ordinary nonlinear least-squares fit. tda_freg stays the more flexible tool for the cases nlreg cannot express at all: any summed objective, not just a sum of squared residuals (a robust loss, a weighted or constrained one, anything written by hand).

Value

An object of class tda_fit. With residuals = TRUE, also carries $residuals.

Writing the function

fn is the quantity summed over cases and minimised, not the regression function itself. Least squares is therefore the squared residual:

tda_freg(c("r = y - a * exp(b * x)", "fn = r*r"), d,
         start = c(a = 1, b = 0.1))

Algorithm

TDA's default minimiser is Newton, which overflows exp() on the first step for functions of this shape. tda_freg therefore defaults to BFGS, which is stable and reaches the same estimates as nls. Pass control = tda_control(algorithm = 5) for TDA's default.

See Also

Other regression: TDA_FAMILIES, TDA_QRMODELS, tda_gdf(), tda_l1reg(), tda_lsreg(), tda_mlrc_design(), tda_mreg(), tda_nlreg(), tda_npreg(), tda_zreg()

Examples

set.seed(1)
d <- data.frame(x = 1:20, y = 3 * exp(0.1 * (1:20)) + rnorm(20, 0, 0.5))
tda_freg(c("r = y - a * exp(b * x)", "fn = r*r"), d,
         start = c(a = 1, b = 0.1))

# the same model, written in R syntax instead
tda_freg({
    r = y - a * exp(b * x)
    fn = r * r
}, d, start = c(a = 1, b = 0.1))

Description

tda_freq counts occurrences of each combination of values across one or more variables; tda_freq1 counts occurrences of each value of one variable; tda_freq2 cross-tabulates two. TDA's manual documents all three (freq, freq1, freq2) on one shared page, since they mostly share the same options (maxcat=/fmt=/tfmt=/df=) and differ mainly in shape. Their R outputs are shaped differently too – tda_freq/tda_freq1 return a table, tda_freq2 a matrix of counts – but all three describe a frequency distribution, so they sit together here as TDA's own manual does.

Usage

tda_freq(..., maxcat = NULL, dir = tempfile("tda"))

tda_freq1(..., maxcat = NULL, options = list(), dir = tempfile("tda"))

tda_freq2(
  ...,
  contingency = FALSE,
  maxcat = NULL,
  options = list(),
  dir = tempfile("tda")
)

Arguments

Value

An object carrying a table for tda_freq/ tda_freq1, a matrix for tda_freq2, and, with contingency = TRUE, a measures named vector (chi-square, Cramer's V, lambda, Somers' D, and the rest).

See Also

Other descriptive statistics: tda_brr(), tda_cov(), tda_dstat(), tda_independence(), tda_ineq(), tda_loglin(), tda_quant(), tda_rcorr(), tda_segr(), tda_subm()

Examples

set.seed(48)
d <- data.frame(g = sample(1:3, 80, TRUE), h = sample(1:2, 80, TRUE))
tda_freq1(d["g"])$table
tda_freq(d)$table  # a joint table over both variables at once, unlike
                   # tda_freq1's single-variable count
tda_freq2(d)$matrix  # a cross-tabulation, rows by g, columns by h
tda_freq2(d, contingency = TRUE)$measures[c("Cramer's V", "Phi")]

Description

frml: like tda_fml, but for episode (survival) data rather than plain cross-sectional rows – fn is the log-likelihood contribution of one episode, and can refer to time (the position within it) the way tda_rate's own define can. This is not fml under another name: the same definitions that work with fml on plain rows fail here with “no episode data defined” until an actual episode structure (survival::Surv()) is given.

Usage

tda_frml(
  formula,
  data,
  definitions,
  start = NULL,
  id = NULL,
  spell = NULL,
  control = NULL,
  constraints = NULL,
  residuals = FALSE,
  residual_vars = NULL,
  protocol = FALSE,
  options = list(),
  dir = tempfile("tda"),
  ...
)

Arguments

Details

fn's expression cannot reference the raw duration/status columns from data directly – Surv()'s columns are consumed by edef() into TDA's internal episode representation, and are not separately visible by name afterward. An fn written against the original column names fails with those names listed as parameters to estimate (TDA's convention for any identifier it does not otherwise recognise), not as a lookup error. Use TDA's built-in accessors for an episode instead: ts and tf (start and finish time), org and des (origin and destination state), time (tf - ts, and the same built-in tda_rate's define can use). A covariate, by contrast, is referenced by its own name from formula, translated automatically the same way tda_fml does it.

Value

An object of class tda_fit. With residuals = TRUE, also carries $residuals.

See Also

Other sequence analysis: tda_evalfi(), tda_fml(), tda_seq_info(), tda_seqgc(), tda_seqm(), tda_seqmd(), tda_seqpe(), tda_seqpm()

Examples

set.seed(1)
n <- 200
d <- data.frame(x = rnorm(n))
haz <- exp(0.5 + 0.8 * d$x)
d$dur <- rexp(n, haz)
d$status <- 1
# the exponential hazard model, written out by hand -- matches
# tda_rate(Surv(dur, status) ~ x, d, model = "exponential"). des and
# (tf - ts) are TDA's built-ins, not d's dur/status columns -- see
# Details.
fit <- tda_frml(Surv(dur, status) ~ x, d,
                c("xb = b0 + x * b1", "gamma = exp(xb)",
                  "fn = des * xb - gamma * (tf - ts)"))
coef(fit)
coef(tda_rate(Surv(dur, status) ~ x, d, model = "exponential"))

# the same model, written in R syntax instead
fit2 <- tda_frml(Surv(dur, status) ~ x, d, {
    xb = b0 + x * b1
    gamma = exp(xb)
    fn = des * xb - gamma * (tf - ts)
})
coef(fit2)

Description

tda_g reaches any of TDA's forty graph commands; the named functions are the ones with an obvious meaning.

Usage

tda_g(g, cmd, options = list(), dir = tempfile("tda"), ...)

Arguments

Value

An object carrying a table and the run.

See Also

Other graph analysis: plot.tda_graph(), tda_dmet(), tda_g_analyses, tda_graph(), tda_locate_line(), tda_ptree(), tda_ucl()

Examples

e <- data.frame(from = c(1, 1, 2, 3), to = c(2, 3, 4, 4), value = 1)
g <- tda_graph(e, directed = FALSE)
tda_g(g, "gmst")$table  # minimum spanning tree

Description

One function per graph command, each a thin wrapper: the command name, the options it takes and the table TDA writes. Everything TDA offers that has no named wrapper is still reachable with tda_g.

Usage

tda_g_components(g, edges = FALSE, ...)

tda_g_shortest(g, format = c("reachable", "matrix", "matrix_no_labels"), ...)

tda_g_mst(g, algorithm = c("kruskal", "prim"), max = FALSE, sort = FALSE, ...)

tda_g_cliques(
  g,
  algorithm = c("bron_kerbosch", "harary_ross"),
  min_size = NULL,
  sort = FALSE,
  ...
)

tda_g_cycles(
  g,
  format = c("fundamental_nodes", "fundamental_edges", "all_v1", "all_v2"),
  ...
)

tda_g_cutpoints(g, per_node = FALSE, ...)

tda_g_flow(g, matrix = FALSE, ...)

tda_g_degrees(g)

tda_g_edges(
  g,
  format = c("edges", "lower_triangle", "lower_triangle_diag", "full", "full_square"),
  ...
)

tda_g_transitive(g, ...)

tda_g_aggregate(g, ...)

tda_g_symmetric(g, ...)

tda_g_toposort(g, relabel_graph = FALSE, ...)

tda_g_nodes(g, ...)

tda_g_dot(g, ...)

tda_g_links(g, backward = FALSE, ...)

tda_g_dcycles(
  g,
  format = c("per_cycle", "per_node", "progress"),
  stop_at_first = FALSE,
  ...
)

tda_g_dblocks(g, ...)

tda_g_compact(g, per_node = FALSE, ...)

tda_g_gcliques(g, ...)

tda_g_random(g, ...)

tda_g_ownership(g, ...)

tda_g_backward(g, ...)

tda_g_paths(g, ...)

tda_g_flowcontrol(g, ...)

tda_g_centred(g, ...)

tda_g_spantrees(g, edges = FALSE, ...)

tda_g_reachable(
  g,
  format = c("count", "list", "strong_nodes", "strong_edges"),
  ...
)

tda_g_neighbourhoods(g, ...)

tda_g_independent(g, min_size = NULL, use_complement = FALSE, ...)

tda_g_spanning(g)

tda_g_eigen(g, n = NULL, tol = NULL, control = NULL, ...)

Arguments