Statalist - Forums

The annoying "mean of xxx" in graph bar

Chen Samulsion — Fri, 06 Nov 2026 10:17:24 GMT

Dear Stata users,

When drawing a graph bar in Stata, if I specify the ascategory option, the label of the variable on the bar chart becomes "mean of xxx." Since I generally don't use graph bar this way—or more often use -statplot- (SSC) as an alternative—I have almost never treated this as a problem. Perhaps because I haven't used graph bar for a long time, when I needed it today for work, I found that this annoying "mean of xxx" cannot be removed through any option setting. Did I overlook some syntax explanation in the help file, or is this simply how it is? Thank you.
P.S. -statplot- can produce very beautiful graph, below is just examples to echo graph bar, and graph combine shares part of the blame.

Code:

sysuse nlsw88
graph bar (mean) wage hours ttl_exp, over(race) ascat name(bar1, replace)
graph bar (mean) wage hours ttl_exp, over(race) name(bar2, replace)
statplot wage hours ttl_exp, over(race) recast(bar) xpose name(bar3, replace)
statplot wage hours ttl_exp, over(race) recast(bar) name(bar4, replace)
graph combine bar1 bar2 bar3 bar4

Array

-ssc new- issues?

Stephen Jenkins — Fri, 06 Nov 2026 09:31:25 GMT

Since yesterday, I've been having problems with -ssc new-

Code:

. ssc new
(contacting http://repec.org)
file http://repec.org/docs/ssc_results.smcl not found
r(601);

Is this related to the issues that arose with php code fixing? Or something else? (I ask because I'm sure -ssc new- worked for me a couple of evenings ago. I used it to download some packages)

linearity assumption check of continuous variables in logistic regression.

seethal davis — Thu, 05 Nov 2026 02:21:55 GMT

hi, i need a help with my logistic regression with binary outcome, where i have to check the linearity assumption of continuous variables with the log-odds of the outcome variable. how should i proceed?

Latest (04 November 2026) update to StataNow 19

Kristin MacDonald (StataCorp) — Wed, 04 Nov 2026 17:52:26 GMT

A new update to StataNow 19 is now available (as of 04 November 2026). This update includes new features for nonlinear heterogeneous DID for panel data, jackknife-after-bootstrap standard errors, forest plot enhancements, find and replace in Word documents, and more. You can read about all of the new features in our latest blog post.

optimize()

Matt Piercy — Wed, 04 Nov 2026 12:10:21 GMT

Hi Stata users,
I am trying to get mata to use a formula for REML from a textbook using the K commutation matrix Array which is different from the process used in mixed. I have got it to work on a simple dataset but when I test it with optimise() I keep getting errors due to an argument mismatch. Do file attached.
Appreciate any help. GitHub couldn't solve it.
Thanks,
Matt

95% CIs from POLYCORR (ssc) do not match 95% CIs from CORRCI (ssc)

Bruce Weaver — Sun, 01 Nov 2026 21:49:46 GMT

The POLYCORR package caught my eye this afternoon as I browsed recent contributions to SSC.

Code:

ssc whatsnew
ssc describe polycorr

It looks potentially useful, but I thought it might be wise to compare its CIs for Pearson correlations to the CIs from Nick Cox's CORRCI package (ssc). And to this point, I have not been able to coax POLYCORR into giving the same results as CORRCI. This strikes me as problematic, so I thought I should alert others. Here is what I have tried so far.

Code:

. // ssc install corrci 
. // ssc install polycorr 
. clear

. sysuse auto 
(1978 automobile data)

. summarize trunk weight length

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
       trunk |         74    13.75676    4.277404          5         23
      weight |         74    3019.459    777.1936       1760       4840
      length |         74    187.9324    22.26634        142        233

. 
. corrci trunk weight length, format(%9.4f)

(obs=74)

               correlations and 95% limits
trunk  weight     0.6722   0.5242   0.7808
trunk  length     0.7266   0.5973   0.8191
weight length     0.9460   0.9154   0.9657

. corrci trunk weight length, format(%9.4f) fisher

(obs=74)

               correlations and 95% limits
trunk  weight     0.6722   0.5209   0.7790
trunk  length     0.7266   0.5940   0.8175
weight length     0.9460   0.9143   0.9653

. 
. polycorr trunk weight length

Polychoric, polyserial and Pearson correlations   Number of obs   =         74
Deletion: listwise                                Variables       =          3
VCE: robust                                       Pairs           =          3
------------------------------------------------------------------------------
Variable 1   Variable 2   Type      rho    Std. err.    z    [95% CI]
------------------------------------------------------------------------------
trunk        weight       Pearson    0.6722  0.0614   10.95   0.5518   0.7926
trunk        length       Pearson    0.7266  0.0593   12.26   0.6104   0.8428
weight       length       Pearson    0.9460  0.0147   64.49   0.9173   0.9748
------------------------------------------------------------------------------
Bivariate normality: no testable polychoric pairs (2x2 tables have df = 0)
Smallest eigenvalue of R =  0.050578  (positive definite)

. polycorr trunk weight length, cifz 

Polychoric, polyserial and Pearson correlations   Number of obs   =         74
Deletion: listwise                                Variables       =          3
VCE: robust                                       Pairs           =          3
------------------------------------------------------------------------------
Variable 1   Variable 2   Type      rho    Std. err.    z    [95% CI, Fisher z]
------------------------------------------------------------------------------
trunk        weight       Pearson    0.6722  0.0614   10.95   0.5336   0.7756
trunk        length       Pearson    0.7266  0.0593   12.26   0.5885   0.8235
weight       length       Pearson    0.9460  0.0147   64.49   0.9085   0.9684
------------------------------------------------------------------------------
Bivariate normality: no testable polychoric pairs (2x2 tables have df = 0)
Smallest eigenvalue of R =  0.050578  (positive definite)

I tried to tag the author of POLYCORR (Yusuf Toyin Yusuf), but he is not among the Yusufs listed, which suggests he is not a Statalist member.

APPENDIX
Here is my code without the output, in case anyone wants it.

Code:

// ssc install corrci 
// ssc install polycorr 
clear
sysuse auto 
summarize trunk weight length

corrci trunk weight length, format(%9.4f)
corrci trunk weight length, format(%9.4f) fisher

polycorr trunk weight length
polycorr trunk weight length, cifz

SSC Archive activity, October 2026

KitBaum — Sun, 01 Nov 2026 14:34:55 GMT

As sec whatsnew shows, a record 49 packages were added to the SSC Archive in October, with multiple contributions from several authors. 39 packages were revised, so running ado update is a good idea. Statitics on package access are available from sec hot, and listings of over 1,600 authors are available from http://repec.org/docs/sscAuthors.html

Sorting Chinese names by surname stroke order in Stata

Chen Samulsion — Sun, 01 Nov 2026 02:37:18 GMT

Dear Stata users,

I need to sort a dataset of Chinese names in Stata according to the standard Chinese surname stroke-order rules (姓氏笔画排序), not just by Unicode or pinyin. I'm hoping someone here has tackled this or can point me to a practical approach.

The sorting rules I need to implement are roughly:

1.Primary key: stroke count of the surname. Fewer strokes come first, more strokes later.

2.Tie-break 1: starting stroke shape (起笔笔形). If stroke counts are equal, compare by the first stroke in the order horizontal (横), vertical (竖), left-falling (撇), dot (点), turning (折). If the first stroke is the same, continue comparing subsequent stroke shapes.

3.Tie-break 2: character structure (字形结构). If stroke count and stroke shapes are both identical, order by left-right structure, then top-bottom structure, then whole-character structure.

4.Same surname: compare given names. If surnames are identical, compare the second character, then the third character of the given name using the same stroke/stroke-order rules.

5.Single-character given names come first. A two-character full name is treated as having a "zero-stroke" second given-name character, so it sorts before three-character full names with the same surname.

6.Compound surnames (复姓): use the first character of the compound surname for sorting. If the first characters are identical, then compare the second character of the compound surname.

So this is more than a simple "count strokes of the surname" sort — it involves stroke order, stroke shape, character structure, given-name tie-breaking, and special handling of single-character given names and compound surnames.

What I've tried / considered:

1.Stata's sort and gsort only sort by numeric values or string collation order (Unicode/pinyin), so they can't directly handle stroke-based ordering.

2.I could pre-compute a numeric "sort key" outside Stata (e.g., in Python with a CJK library) and import it, but I'd prefer a reproducible workflow if possible.

My questions:

1.Is there any Stata package, ssc module, or user-written command that supports Chinese stroke-order sorting?

2.If not, what's the most practical reproducible workflow? For example, generating a full sort key in Python/R and importing it into Stata?

3.Has anyone implemented the full stroke-order + stroke-shape + structure rules (e.g., via Unihan data, cjklib, or similar) and then used the resulting rank as a sort variable in Stata?

4.Any pitfalls with compound surnames, single-character given names, or rare characters I should watch out for?

I'm working with lists of about 100 names each, so performance isn't a concern — correctness and reproducibility are. Any pointers, code snippets, or package recommendations would be greatly appreciated.

Thanks!

Code:

* Example generated by -dataex-. To install: ssc install dataex
clear
input str9 name double(raworder strokeorder)
"路甬祥"  1 75
"王佛松"  2  4
"姚平录"  3 50
"汪前进"  4 27
"吴文俊"  5 23
"潘吉星"  6 78
"王渝生"  7  7
"韩琦"     8 73
"杜石然"  9 17
"薄树人" 10 79
"陆敬严" 11 32
"王兆春" 12  3
"梅汝莉" 13 68
"金秋鹏" 14 37
"徐苹芳" 15 56
"谈德颜" 16 66
"师昌绪" 17 11
"俞伟超" 18 48
"郭湖生" 19 59
"戴念祖" 20 80
"罗桂环" 21 36
"赵承泽" 22 42
"马素卿" 23  2
"梁家勉" 24 70
"郭书春" 25 58
"黄汲清" 26 67
"艾素珍" 27  8
"钱临照" 28 55
"吴瑰琦" 29 24
"谭其骧" 30 77
"王大珩" 31  1
"钱三强" 32 52
"王绶琯" 33  6
"唐寰澄" 34 64
"孙枢"    35 14
"严东生" 36 16
"曾世英" 37 74
"杨熺"    38 22
"陈美东" 39 33
"宋正海" 40 29
"丘光明" 41  9
"蒋顺学" 42 72
"刘钝"    43 13
"钱伟长" 44 54
"章综"    45 69
"廖育群" 46 76
"唐锡仁" 47 63
"李家治" 48 20
"席泽宗" 49 61
"李家明" 50 19
"侯仁之" 51 47
"张存浩" 52 30
"席龙飞" 53 60
"周嘉华" 54 40
"汪子春" 55 26
"涂光炽" 56 65
"柯俊"    57 45
"汪德昭" 58 28
"周光召" 59 38
"杜石然" 60 18
"范楚玉" 61 35
"钱文藻" 62 53
"徐冠仁" 63 57
"姜丽蓉" 64 49
"孙鸿烈" 65 15
"张含英" 66 31
"王振铎" 67  5
"袁翰青" 68 51
"白寿彝" 69 10
"武衡"    70 34
"胡启恒" 71 43
"胡道静" 72 44
"华觉明" 73 12
"董恺忱" 74 71
"赵匡华" 75 41
"席泽宗" 76 62
"柯俊"    77 46
"杨文衡" 78 21
"余志华" 79 25
"周魁一" 80 39
end

label var name "Chinese name"
label var raworder "Raw order of names"
label var strokeorder "Stroke order of names"

replacevalue of variable

Euslaner — Sat, 31 Oct 2026 17:19:45 GMT

I have data as follows:

fsum Heritage_Reformed_Churches International_Churches_of_Christ Moravian_Church_in_America__Sout National_Association_of_Congrega - Open_Bible_Standard_Churches

Variable | N Mean SD Min Max
----------------------------------+---------------------------------------------
Heritage_Reformed_Churches | 409 1.97 22.59 0.00 412.00
International_Churches_of_Christ | 409 97.51 126.94 0.00 752.00
Moravian_Church_in_America__Sout | 409 9.26 53.60 0.00 573.00
National_Association_of_Congrega | 409 158.95 350.52 0.00 2627.00
National_Association_of_Free_Wil | 409 506.29 1417.55 0.00 11653.00
National_Baptist_Convention_of_A | 409 706.51 1353.97 0.00 9978.00
National_Baptist_Convention__USA | 409 4394.00 7272.73 0.00 58446.00
National_Missionary_Baptist_Conv | 409 608.95 906.11 0.00 5385.00
Non_denominational | 409 27897.60 13623.09 1864.00 71184.00
North_American_Baptist_Conferenc | 409 135.31 503.92 0.00 5220.00
Old_Order_River_Brethren | 409 1.51 19.94 0.00 360.00
Open_Bible_Standard_Churches | 409 53.37 167.78 0.00 1330.00

I typed:

replace evangelicalnewrev = 1 if Heritage_Reformed_Churches > 0

evangelicalnewrev has been define earlier. However I received: (0 real changes made)

There must be a ay to change the value of evangelicalnewrev since Heritage_Reformed_Churches has most of it values > 0, How can I do this?

.

New package on SSC: tvtools - build time-varying exposure data for survival analysis

Tim Copeland — Fri, 30 Oct 2026 23:53:57 GMT

Thanks to Kit Baum, a new package, tvtools, is now available on SSC. It requires Stata 16 or later.

tvtools turns prescription, treatment and other episode records into analysis-ready time-varying survival data. Anyone who has done this by hand in a pharmacoepidemiology study knows the work: overlapping prescriptions, gaps that should or should not count as exposure, follow-up that must be tiled without losing or double-counting a day, several drugs merged onto one timeline, and outcomes and competing events that must cut the intervals at the right place. tvtools does each of these steps as a command, with checks at each step. I was inspired to make this command because after doing everything manually, I never wanted to do it again.

The commands:
- tvexpose: episode records to person-time intervals, with categorical, ever/never, current/former, cumulative-duration, recency and cumulative-dose definitions, plus grace periods, lags and washout windows
- tvmerge: combine several time-varying exposures into one set of intervals
- tvevent: attach single, competing or recurrent events, with gap-time clocks for recurrent events
- tvage, tvband, tvsplit: split follow-up by exact age, calendar time, elapsed time, or several Lexis timescales at once
- tvpanel: entry-anchored person-period panels, for example for marginal structural models
- tvbuild and tvspec: the whole construction as one validated call, with a provenance record
- tvdiagnose: coverage, gaps, overlaps, exposure summaries and swimlane plots
- tvweight: inverse probability of treatment weights (optionally stabilized), overlap and matching weights, and inverse probability of censoring weights, with balance and positivity diagnostics

Dialog boxes are included for tvexpose, tvmerge and tvevent.

A small example: four people, prescription episodes for two drug classes, and an outcome with death as a competing event.

Code:

* cohort: one row per person
clear
input long id study_entry study_exit byte female
1 21915 22280 1
2 21915 22280 0
3 21930 22280 1
4 21930 22280 0
end
format study_entry study_exit %td
tempfile cohort episodes events
save `cohort'

* prescription episodes
clear
input long id rx_start rx_stop byte rx_class
1 21919 21965 1
1 22024 22069 2
2 22000 22040 1
3 21950 22030 2
end
format rx_start rx_stop %td
label define rx 0 "Unexposed" 1 "SSRI" 2 "SNRI"
label values rx_class rx
save `episodes'

* outcome and death dates
clear
input long id event_date death_date
1 22010     .
2     .  22100
3 22020     .
4     .     .
end
format event_date death_date %td
save `events'

* 1. episodes to exposure intervals covering all of each person's follow-up
use `cohort', clear
tvexpose using `episodes', id(id) start(rx_start) stop(rx_stop) exposure(rx_class) reference(0) entry(study_entry) exit(study_exit) frameout(f_drug)

* 2. cut the intervals at the outcome, with death as a competing event
use `events', clear
tvevent, frame(f_drug) id(id) date(event_date) compete(death_date) start(rx_start) stop(rx_stop) generate(status)
list id rx_start rx_stop tv_rx_class status, sepby(id) noobs

* 3. stset: tvtools intervals are closed [start, stop] dates, so time0 is start - 1
generate double t0 = rx_start - 1
stset rx_stop, id(id) failure(status == 1) time0(t0)

The listing shows each person's follow-up tiled by unexposed, SSRI and SNRI periods, ending at the outcome (status 1, labelled "Event"), at death (status 2, "Competing"), or at the end of follow-up (status 0, "Censored"). From there it is ordinary stcox or stcrreg, or my finegray package.

Some design points:
1. Every day of follow-up is accounted for. Intervals are closed, inclusive daily dates that tile each person's follow-up, with the reference category filling uncovered time.
2. Quantities add up. Dose, rate and cumulative variables are apportioned by days when intervals are split, so totals still sum to the source.
3. Identifiers are never guessed: source and cohort IDs must agree in name and type, and the commands never pad, recast or remap them. Person-time lost when merging exposures that do not cover the same days is reported, not dropped silently.
4. tvbuild runs a whole construction (one or several sources, with events) as one call: it checks the plan before changing anything, builds through the same engines as the individual commands, and commits the result frame and a provenance manifest together. Its dryrun option checks the plan without building anything.

To install:

Code:

ssc install tvtools
help tvtools

The source, a README with worked examples of every command, and an issue tracker are at https://github.com/tpcopeland/Stata-...e/main/tvtools. Bug reports and suggestions are welcome, here or on GitHub.

New package on SSC: rangematch - range and interval-overlap joins using frames

Tim Copeland — Fri, 30 Oct 2026 23:18:30 GMT

Thanks to Kit Baum, a new package, rangematch, is now available on SSC. It requires Stata 16.1 or later.

rangematch joins the data in memory to a dataset on disk or in a frame by ranges rather than exact keys: each master observation is paired with every using observation whose key falls inside the master's [low, high] interval, or, with overlap(), every using interval that overlaps it. Typical uses are events inside exposure windows, lab values within a period around a visit, or prescriptions overlapping follow-up. The alternative, joinby followed by a keep if, builds every within-group pair first and can exhaust memory on registry-sized data.

Adverse events that fall inside each patient's drug-exposure windows:

Code:

* exposure windows: one row per patient and drug
clear
input long patient_id str10 start_string double exposure_days str9 drug
101 "2020-01-15" 30 "drug_a"
101 "2020-03-01" 14 "drug_b"
102 "2020-02-10" 21 "drug_a"
103 "2020-02-20" 10 "drug_c"
end
generate double exposure_start = daily(start_string, "YMD")
generate double exposure_end = exposure_start + exposure_days
format exposure_start exposure_end %td
drop start_string exposure_days
tempfile exposures events
save `exposures'

* adverse events: one row per event
clear
input long patient_id double event_id str10 event_string str12 event_type
101 1001 "2020-01-20" "rash"
101 1002 "2020-02-20" "headache"
101 1003 "2020-03-10" "nausea"
102 1004 "2020-02-15" "dizziness"
102 1005 "2020-03-20" "fatigue"
103 1006 "2020-02-25" "rash"
end
generate double event_date = daily(event_string, "YMD")
format event_date %td
drop event_string
save `events'

* events inside each exposure window, same patient only, written to a frame
use `exposures', clear
rangematch event_date exposure_start exposure_end using `events', by(patient_id) keepusing(event_id event_date event_type) frame(exposure_events) stats
frame exposure_events: list patient_id drug exposure_start exposure_end event_date event_type, sepby(patient_id) noobs

Each window picks up the one event inside it; events 1002 and 1005 fall outside every window and are not joined. The exposure data in memory are unchanged.

The main features:
1. Bounds are variables, scalar offsets from the key (rangematch event_date -30 30 using ...), or . for an open bound. closed() sets which endpoints are included, tolerance() absorbs floating-point noise at the boundaries, and missing() decides whether a missing bound is open, dropped or an error.
2. overlap(ulow uhigh) matches master intervals to every overlapping using interval, for example follow-up periods against prescription episodes.
3. nearest() with ties() keeps only the closest match before, after, or in either direction, and masterid(), usingid() and distance() record where each output row came from and how far apart the keys are.
4. Output replaces the data in memory, or goes to a frame() or a dataset on disk with saving(). unmatched() controls whether unmatched master and using rows are kept, and generate() adds a merge-style indicator.
5. dryrun reports how many pairs a join would produce without producing them, and maxpairs() stops before an oversized join is materialized. assert() checks that required matches exist.
6. Pairs are generated by binary-search and sweep algorithms in Mata, so the Cartesian product is never built. There are no external dependencies.

How rangematch differs from rangejoin
rangematch follows the bound syntax of Robert Picard's rangejoin, and by(), keepusing(), prefix(), suffix() and all work the same way, so a call such as rangejoin key lo hi using file, by(group) runs unchanged as rangematch key lo hi using file, by(group). Compared with rangejoin 1.1.3 (SSC, April 2021):

Code:

                       rangejoin 1.1.3            rangematch
---------------------------------------------------------------------------------
Stata version          11 or later                16.1 or later
Dependencies           rangestat (SSC)            none
Key inside interval    yes                        yes
Interval overlap       no                         overlap(ulow uhigh)
Endpoints              always [lo, hi]            closed(): [lo, hi], [lo, hi), (lo, hi], (lo, hi)
Boundary tolerance     no                         tolerance()
Missing bound          open-ended                 missing(): open-ended, drop, or error
Nearest match only     no                         nearest(), ties()
Provenance             no                         distance(), masterid(), usingid()
Using data             dataset on disk            dataset on disk or a frame
Output                 replaces data in memory    data in memory, frame(), or saving()
Unmatched rows         master rows always kept    unmatched(master|using|both|none)
Match indicator        no                         generate()
Size check before      no                         dryrun, maxpairs()
Required matches       no                         assert()
---------------------------------------------------------------------------------

Two rangejoin behaviours do not carry over. First, rangejoin key 0 10 using file treats 0 and 10 as absolute bounds when the master data have no variable named key; rangematch reads scalar bounds only as offsets from the master key and stops with r(111), so generate lo and hi first. Second, rangejoin never matches a master row whose key is missing, whereas rangematch with variable bounds matches on lo and hi alone, so such a row can match. The help file ("Migrating from rangejoin") gives the recipe for both.

With about one match per master row the two are close, and rangejoin is faster on the smallest join; rangematch pulls ahead as the number of pairs grows. Times in seconds on synthetic point-in-interval joins with equal master and using sizes (Stata/MP 17 on Linux, single run; both commands returned the same pairs; rangematch run with unmatched(none) nosort):

Code:

rows        matches per row        pairs   rangematch   rangejoin   rangejoin/rangematch
-----------------------------------------------------------------------------------------
10,000      1                     10,000        0.151       0.084                   0.56
10,000      up to 21             207,800        0.216       0.345                   1.60
100,000     1                    100,000        0.540       0.508                   0.94
100,000     up to 11           1,098,500        1.053       2.365                   2.25
1,000,000   1                  1,000,000        4.005       4.415                   1.10
1,000,000   up to 3            2,999,800        5.263      10.012                   1.90
-----------------------------------------------------------------------------------------

The benchmark script, bench_rangematch.do, is in the GitHub repository linked below if you want to rerun it on your own machine.

If you need summary statistics over a range rather than the joined rows themselves, rangestat (Picard, Cox and Ferrer) is the better tool.

To install:

Code:

ssc install rangematch
help rangematch

The source, a README with worked examples, and an issue tracker are at https://github.com/tpcopeland/Stata-...ain/rangematch. Bug reports and suggestions are welcome, here or on GitHub.

New package on SSC: kmplot - Kaplan-Meier and cumulative incidence curves with risk tables

Tim Copeland — Fri, 30 Oct 2026 22:09:47 GMT

Thanks to Kit Baum, a new package, kmplot, is now available on SSC. It requires Stata 16 or later.

kmplot draws Kaplan-Meier survival or cumulative failure curves from stset data, or Aalen-Johansen cumulative incidence curves when there is a competing event, with confidence bands, a number-at-risk table, median lines, censoring marks and a log-rank p-value, in one command.

Code:

sysuse cancer, clear
label define type 2 "Drug A" 3 "Drug B", modify
stset studytime, failure(died)
kmplot, by(drug) ci median medianannotate pvalue censor risktable riskevents timepoints(0 5 10 15 20 25 30 35)

(The label define line only names the two active drugs, which cancer.dta labels "Other" and "NA".)

With a competing event, compete() gives the cumulative incidence of the stset failure, with the competing event taken out of the risk set rather than censored:

Code:

* time to SI phenotype after HIV infection, with AIDS as the competing event
webuse hiv_si, clear
stset time, id(patnr) failure(status == 2)
kmplot, by(ccr5) compete(status == 1) ci risktable tmax(10) xtitle("Years from HIV infection")

The main features:
1. The risk table shares the plot's time axis, and its counts respect delayed entry and multiple-record data. riskevents adds cumulative events, and riskmask(5) prints counts below 5 as <5 for registries with disclosure rules. With probability weights the table counts people by default; riskweight() changes that.
2. Confidence intervals follow sts conventions: Greenwood standard errors with log-log (the default), log or plain intervals, drawn as shaded bands or as lines. For cumulative incidence without weights, the standard errors and ln(-ln) bounds match stcompet.
3. median adds median survival lines, landmark() reports estimates at fixed times with confidence limits, and pvalue adds the log-rank test. The numbers are returned in r(medians), r(landmarks) and r(p), so the figure and the text of the paper use the same values.
4. Curves use the Okabe-Ito palette by default, which stays distinguishable for readers with red-green colour-vision deficiency for up to five groups. colors() takes Stata colour names or RGB triplets, and ordinary twoway options pass through.
5. failure plots 1 - S(t); percent, tmax() and byorder() set a percent axis, the plotted horizon and the group order. saving() and risksaving() write the curve and risk-table data to datasets, and export() writes the graph.

All of this can be assembled from sts graph, sts list and stcompet; kmplot puts it in one call with publication defaults and returns the numbers behind the picture. For the regression side of competing risks, see my finegray package, also on SSC.

To install:

Code:

ssc install kmplot
help kmplot

The source, a README with a gallery of examples, and an issue tracker are at https://github.com/tpcopeland/Stata-...ee/main/kmplot. Bug reports and suggestions are welcome, here or on GitHub.

Controlling group order when using egen group()

Alex Greene — Fri, 30 Oct 2026 22:03:44 GMT

I am creating a single numeric group ID from three string variables with egen group(). That works, but the new IDs are numbered in alphabetical order of the combinations. I would like group 1 to be the most frequent combination, group 2 the next most frequent, and so on. Right now I get the frequencies with collapse, sort the result descending and merge it back by hand. Is there an easier way to control that ordering, for example an option in egen group() or a built-in command? I am using Stata 18 on Windows. Thanks for any suggestions.

New package on SSC: eplot - forest plots and coefficient plots from data, estimates, matrices or frames

Tim Copeland — Fri, 30 Oct 2026 21:58:11 GMT

Thanks to Kit Baum, a new package, eplot, is now available on SSC. It requires Stata 16 or later.

eplot draws forest plots and coefficient plots with one command, whatever form the estimates are in:
- variables in memory, such as study-level results for a meta-analysis or subgroup hazard ratios typed in from a paper
- the active estimation results, or several stored models side by side
- a matrix of estimates and standard errors, or estimates and confidence limits
- a graph-ready frame

eplot works out which of these you mean from the call, and the options for selecting, labelling, transforming and displaying the estimates are the same in every case.

A coefficient plot of the model you just fitted, and two stored models compared:

Code:

sysuse auto, clear
logit foreign mpg weight length headroom
eplot ., drop(_cons) eform values

quietly regress price mpg weight foreign
estimates store base
quietly regress price mpg weight length headroom foreign
estimates store extended
eplot base extended, drop(_cons) modellabels("Base" "Extended") cicap

A meta-analysis layout from results held as data, with weighted boxes, prediction intervals, a column of formatted estimates, a pooled diamond, heterogeneity text and direction labels:

Code:

clear
input str20 study double(es lci uci pi_lci pi_uci weight) byte type
"Smith 2018"   -0.42  -0.78  -0.06  -1.15   0.31  12.3  1
"Jones 2019"   -0.31  -0.58  -0.04  -1.04   0.42  16.8  1
"Brown 2020"   -0.18  -0.41   0.05  -0.91   0.55  21.5  1
"Lee 2021"     -0.55  -0.93  -0.17  -1.28   0.18  10.2  1
"Garcia 2022"  -0.27  -0.49  -0.05  -1.00   0.46  19.1  1
"Patel 2023"   -0.09  -0.35   0.17  -0.82   0.64  20.1  1
"Overall"      -0.28  -0.41  -0.15   .       .      .    5
end

eplot es lci uci, labels(study) weights(weight) type(type) values vformat(%4.2f) pi(pi_lci pi_uci) i2("42.1%") tau2("0.021") qstat("8.63, df=5, p=0.125") effect("Mean Difference (95% CI)") favors("Favors Treatment" "Favors Control")

eplot draws meta-analysis results; it does not compute them. The pooled estimate and the heterogeneity statistics come from meta, metan, or wherever you estimated them.

Some other points:
1. Confidence intervals match Stata's own tables: eplot uses the t distribution with e(df_r) where the model has it, parameter-specific degrees of freedom after mi estimate and mixed, dfmethod(), and eform exponentiates only the equations that Stata's own eform option would.
2. keep(), drop(), rename(), groups() and headers() organize the rows; eform, logscale, cicap, stars, a values column with its own format and header, labelled reference lines, and favors() labels cover most publication needs. style() applies presets modelled on Lancet, JAMA, NEJM and BMJ figures, and anything you specify overrides the preset.
3. No silent drops: a keep() or drop() name that matches nothing is an error, and a data row with a missing estimate is left out with a note naming it.
4. The plotted estimates, limits and p-values are returned in r(), together with the twoway command that was executed, so the graph can be taken further by hand.
5. With my tabtools package (also on SSC), regtab and effecttab can write a graph-ready frame alongside the table through eplotframe(), so the same model results give both the table and the figure.

If you already use Ben Jann's coefplot, it remains the more comprehensive tool for coefficient plots. eplot's niche is forest plots and coefficient plots with one syntax, from data in memory and frames as well as estimation results.

To install:

Code:

ssc install eplot
help eplot

The source, a README with a gallery of examples, and an issue tracker are at https://github.com/tpcopeland/Stata-...ree/main/eplot. Bug reports and suggestions are welcome, here or on GitHub.

New package on SSC: fvgen - turn factor-variable interactions into labeled variables

Tim Copeland — Fri, 30 Oct 2026 21:42:32 GMT

Thanks to Kit Baum, a new package, fvgen, is now available on SSC. It requires Stata 16 or later.

Factor-variable notation is the right way to fit a model with interactions, but the results are awkward to export: the coefficients are named things like 1.foreign#c.mpg, base levels appear as their own rows, and nothing carries a readable label. fvgen expands a factor-variable specification into ordinary variables, an indicator for each non-base level and a product for each interaction term, labelled from your value labels, so a table-making command gets one clean, labelled row per coefficient.

Code:

sysuse auto, clear
regress price i.foreign##c.mpg
fvgen i.foreign##c.mpg
local flat `r(allvars)'
describe `flat'
regress price `flat'

The second regression is the same model: in the default (uncentered) setup a fit on r(allvars) reproduces the factor-variable model's fit, and when the design has full rank every coefficient and standard error matches as well. The generated variables are labelled Foreign, Mileage (mpg), and Foreign × Mileage (mpg). Interaction labels join their parts with × by default; xsymbol() changes the symbol.

Other options:

Code:

* a different reference level, named in the main-effect labels
sysuse auto, clear
fvgen i.foreign##i.rep78, ref(rep78 3) vsref("(vs. @)")
regress price `r(allvars)'
fvgen, drop

* a separate mpg slope within each level of foreign, instead of a reference slope plus a difference
fvgen i.foreign##c.mpg, simple(foreign)
regress price `r(allvars)'
fvgen, drop

* margins still works: refit with native factor variables and store the clone
fvgen i.foreign##c.mpg
regress price `r(allvars)'
fvgen, margins store(m_price)
estimates restore m_price
margins, dydx(mpg) at(foreign=(0 1))

Also:
1. center mean-centres continuous terms before forming products, alllevels materializes every level including the base, and prefix() sets the name prefix.
2. Nothing is overwritten by accident. Before changing the data, fvgen refuses any plan in which two terms would map to the same generated name or an output would overwrite a source variable.
3. Generated variables carry provenance characteristics, so fvgen, drop removes only the variables fvgen created and leaves your own variables, renamed copies and clonevar copies alone.
4. fvgen, margins refuses (r(498)) if the data have changed since fvgen ran, rather than handing margins a model that no longer matches.

Stata's xi also creates indicators and products. fvgen adds labels built from value labels, modern factor-variable syntax (ib#., ibn., c., ##), per-factor reference levels and per-group slopes, a matching teardown, and the bridge back to margins.

To install:

Code:

ssc install fvgen
help fvgen

The source, a README with worked examples, and an issue tracker are at https://github.com/tpcopeland/Stata-...ree/main/fvgen. Bug reports and suggestions are welcome, here or on GitHub.