Despite an increasing interest in digital subthreshold circuits little research has been dedicated to timing modeling in this voltage domain so far. Especially high timing variabilities makes proper modeling necessary to allow for the prediction of timing behavior and timing yield on the path towards design automation. This paper first deals with gate timing characterization at sub-threshold voltages and a characterization waveform well resembling the actual transistor-level waveforms in this voltage domain is proposed. The error made in this abstraction step is identified and shown to be typically below 3%. Secondly, the modeling of timing variability is considered and the high correlation between gate delays due to slope propagation combined with strong non-linearities in the delay-slope dependencies are pointed out as modeling challenges. A path-based logic-level Monte-Carlo technique, already magnitudes faster than transistor-level simulation, is applied and shown to match transistor-level Monte-Carlo simulation results better than 3% in mean and 7% in standard deviation values.