Sunday, September 13, 2015

Trust me, I am a developer!

So, that is the end. After 5 years with Blogger, I've decided to move onto Ghost and Digital Ocean. In other words - let me introduce trustmeiamadeveloper.com!

I've copied all the English articles from this blog into the new one. Plus, seems finally I have a lot of new topics to share! I do believe that in the next few months you'll find something interesting on my new blog. Trust me, I am a developer! :)
Read more...

Sunday, June 7, 2015

Logistic regression in Excel

Yes, that is weird :) If you need to deal with statistics you have to use a special software like Mathlab or Statistica. But if you are limited in your choice and Excel is the only instrument you have, this manual is for you :)


Linear regression


First of all, Excel already has the "Regression" add-in which allows you to perform a simple lineral regression analysis:

Unfortunately, logistic regression isn't supported by that add-in. But there is always a way to workaround a problem!


Solver


We're starting our journey from an another add-in which name is "Solver". This add-in allows us to solve different minimization/maximization tasks.

Here is an example of a simple maximization problem: we have a furniture factory which produces 2 models of cabinet (A and B). Each model requires a different quantity of resources (wood, time) and generates different income. Also the factory has a limited quantity of available woods and (of course) time. Using simplex-like methods Solver allows us to maximise a target function (the sum of potential income, in current case) by determining the quantity of products for each model that should be produced.

The same idea (solving a maximization task) could be used to compute logistic regression. Let's see how it can be done.


Logistic regression


Imagine that you are an owner of a company and you have a database of clients. Some of them terminated the contract with your company during the last year and you want to predict which clients are thinking about leaving you right now :)

Let's say that you think that decision to stay with your company depends on client's sex and age. Then you can extract the follwing data for the previous year (Trainig Set):

Our next step makes a proposal about how the objective function should look like. In common case it looks like . We assume that there is a linear dependency between decision to leave and sex/age of a client: . In other words, let's just put (theta is a canonical letter for this case, but using it in Excel is a little bit difficult).

The calculation below requires some "initial" values for values. For the moment, let's put them eqauls 1 (the values will differ from the ones you see here).
After that we have to calculate the logit function ;
The main advantage of this function is tending to 0 for x < 0 and tending to 1 for x > 0. So, we can say that and . Taking into account that y belongs to {0, 1}, (Bernoulli distribution).

Hence, our task has been reduced to the selection of the theta () parameters of Z of the objective function Z to maximize P{y|x} probability.

There are two important moments:
  • the Solver can search for a local maximum only, so you have to guess "valid" initial values of variables;
  • the boundary conditions (Y=1,X=1 and X=0,Y=0) must be treated separately (=IF(OR(AND(A5=1,L5=1),AND(A5=0,L5=0)),1,(L5^A5)*(1-L5)^(1-A5)));
Finally, when we found values, we can resote the objective function Z = -25 + 0.586549634 * Age - 1.66748138716445 * Gender.

Let's check the result on our training set:

Where the "Has Terminated" column above is the logit function (ROUND(1/(1+EXP(-T5)),0)).

Now we're sure that our parameters are correct and we can use them to "predict" the future. So, let's do that!



Conclusion


As you can see, the way was a little bit tricky. But anyway, we've managed to pass it! :)


Read more...

Saturday, December 20, 2014

ELK as a software monitoring tool

Development is only the first step of an application lifecircle. Once it is in production (or even in testing) you have to monitor it: watch if errors occur, memory leaks and how fast its components work. Hopefully, there are a lot of solutions which allow you to do it. Today, I want briefly describe my experience with one of them - ELK.


ELK is a set of 3 products created by Elasticsearch company:
  • Elasticsearch - Lucene based distributed search engine;
  • Logstash - utility for data (logs) extraction and parsing;
  • Kibana - visualization tool;
All these components are well described on the related pages, so I don't see any reason to repeat anything here. Instead, let's consider how this technology can be used for your application.



Architecture


Imagine that we have a distributed application running on N hosts. Each host contains a few components of the application which write logs on the local file system (gathering logs on NFS-like system is not a good idea, yeah?).

How can we organize the monitoring process in this case? Well... it's very simple! First of all we need to find a separate host and install Elasticsearch there. Yes, we know that it can be used in the distributed mode... but for the first time let's stop on the simplest solution. On the same host we can install Kibana, since we know that it's a quite lightweight html5 application almost running on the client side.

The final step is installing Logstash on all N host and configuring it to parse certain logs and send them into the Elasticsearch instance (through its REST api).





Data for monitoring


The first question we should answer is what kind of data/metrics we want to monitor? I think we can distinguish at least three areas here:
  • application errors;
  • application metrics: amount of processed tasks, time that spent for major operations, etc;
  • system metrics: cpu usage, heap size, gc pauses (through MBeans and Sigar / a simple external bash script);

The first option seems to be obvious - if your application says it has an exception/error, you have to fix it. We can easily grep such kind of information from the logs using a simple pattern "[ERROR]" (or something like that). But strictly speaking, we also need the date when the error occurred... and it's a good idea to extract the error's description from the message. Hopefully, Logstash knows how to deal with Log4j logs. But what if we want to deal with more complicated, structured data as our metrics are?

For this case there is an extremely old but very effective idea - "if you want your logs to be used not only by humans, write them in some parsable format: CSV, XML or JSON". Sounds great, doesn't it? :)
Of course, Logstash can deal with all these formats. So, all we need is a separate appender for your metrics class with a JSON layout (see below).

After that you'll be able to have all your logs in JSON format! But from my point of view, it better to proceed storing "classic" logs (including errors) in the non-formatted way (which is more acceptable for humans) and use JSON only for metrics (which are mostly supposed for machine analysis).
{
  "@version":1,
  "@timestamp":"2014-01-27T19:52:35.334Z",
  "host":"app1",
  "component":"WEB",
  "type": "SYS",
  "metric": "CPU_USAGE",
  "value": 0.85
}
{
  "@version":1,
  "@timestamp":"2014-01-27T19:52:35.738Z",
  "host":"app1",
  "component":"BACKEND",
  "thread_name":"processor-14",
  "type": "APP",
  "metric": "BACK-SCENARIO1",
  "value": 13212
}
Hint: in a real log each json object must be written in one line.

Taking into account that we have Logstash installed on each node and each component has its own log file, some fields in example above are exceed and can be moved into the Logstash configuration (e.g. "host"). Also keep in mind, that ElasticSearch builds indexes for all the fields you have... so, if you have problems with free space, do not add the data you aren't going to use (or try to amend ElasticSearch behavior).



Installation


You can download tarballs of all the components by the links I put above, but I would recommend you to use the package based way. In this case you wont have to create all these start up scripts, service wrappers and other stuff.



Configuration


Since we aren't going to build Elasticseach cluster, its configuration will be quite simple:
  • update /etc/elasticseacrh/elasticsearch.conf to use actual IP address;
  • configure indexes clean up by adding Curator into CronTab;


Initial Kibana configuration isn't required at all (taking into account the we installed it on the same host as the Elasticsearch server). But if you want to get some profit from using ELK you have to spent a few hours on configuring histograms for your metrics (see "Results" section).

So, the most interesting part is Logstash. It has a very large range of abilities and settings. Let's consider the configuration (/etc/logstash/conf.d/01-your-app.conf), which works with both "classic" and "formatted" logs:
input {
        file {
                start_position => beginning
                path => [ "/var/log/your_app/*.log"]
                sincedb_path => "/var/data/logstash/sincedb"
                type => "file_log4j"
        }
        file {
                start_position => beginning
                path => [ "/var/log/your_app/*.log.json"]
                sincedb_path => "/var/data/logstash/sincedb"
                codec => json
                type => "file_json"
        }
}

filter {
        if [type] == "file_log4j" {

                # Parse Log4j entries
                grok {
                        match => ["message", "%{TIMESTAMP_ISO8601:timestamp}%{SPACE}\[%{DATA:thread}\]%{SPACE}%{LOGLEVEL:level}%{SPACE}%{DATA:classname}%{SPACE}%{GREEDYDATA:message}"]
                        overwrite => [ "message" ]
                }
                
                # Filter out non-ERRORs
                if [level] != "ERROR" {
                        drop {}
                }

                # Update entry's date
                date {
                        match => [ "timestamp", "yyyy-MM-dd HH:mm:ss,SSS" ]
                        remove_field => [ "timestamp" ]
                }
        }
}


output {
        elasticsearch {
                host => elasticsearch.host
        }
}
Now, when our system is configured, it's time to take a look at the results!



Results


Here is an example of a real monitoring system based on ELK. It's too big to paste its whole screenshot into this article, so I'll consider it by parts. Also I slightly GIMPed the photos... just... in case :)

First of all there is a set of filters which are used in diagrams below. Under ordinal circumstances they are folded. Also you always can specify some extra filtering criteria (e.g. time interval).

Below you can find a table with a list of the latest errors. Keep in mind that the data is updated in real time (according to predefined time interval - 5s).

The last stop is a set of histograms which show us different metrics for different components (how many tasks in the queue, how long they were stored in the database, amount of user requests in the web interface, etc).



Going further


As you can see from the text above, storing metrics onto the file system is an exceed step - we don't need them there. Instead we can directly send them into Logstash. This fact leads us to the new architecture:



Here, we can have only one Logstash instance (on the same host as ElasticSearch) listening to UDP socket. In this case we avoid writing data on the local FS and configure our Log4j appender to write directly into Logstash. This approach slightly reminds how SysLog works.

Log4J config:
# Root logger option
log4j.rootLogger=INFO, logstash
 
# Direct log messages to Logstash
log4j.appender.logstash=org.apache.log4j.receivers.net.UDPAppender
log4j.appender.logstash.remoteHost=logstash.host
log4j.appender.logstash.port=6666
log4j.appender.logstash.application=YourProject
log4j.appender.logstash.encoding=UTF-8
log4j.appender.logstash.layout=name.krestjaninoff.log.json.LogstashJsonLayout

Logstash config:
input {
        udp {
                port => 6666
                codec => json
                type => "udp_json"
        }
}

...

output {
        elasticsearch {
                host => localhost
        }
}

LogstashJsonLayout code (just a proof of concept):
/**
 * JSON layout for Logstash
 */
public class LogstashJsonLayout extends Layout {
    
    private final Gson gson = new GsonBuilder().create();
    private final String hostname = getHostname().toLowerCase();
    private final String username = System.getProperty("user.name").toLowerCase();
    
    @Override
    public String format(LoggingEvent le) {
        Map<String, Object> message = new LinkedHashMap<>();
        
        message.put("@timestamp", le.timeStamp);
        message.put("hostname", hostname);
        message.put("username", username);
        message.put("level", le.getLevel().toString());
        message.put("thread", le.getThreadName());
            
        if (le.getMessage() instanceof Map) {
            message.putAll(le.getMessage());
            
        } else {
            message.put("message", null != le.getMessage() ? le.getMessage().toString() : null);
        }
        
        return gson.toJson(message) + "\n";
    }
   
    private static String getHostname() {
        String hostname;
        try {
            hostname = java.net.InetAddress.getLocalHost().getHostName();
        } catch (Exception e) {
            hostname = "Unknown, " + e.getMessage();
        }
        return hostname;
    }

    @Override
    public boolean ignoresThrowable() {
        return false;
    }
    
    @Override
    public void activateOptions() {
    }
}

Read more...

Sunday, November 23, 2014

Using Hidden Markov Models for speech recognition

In the previous article we considered the simplest variant of a speech recognition application. Now I want to suggest to plunge into this topic a little deeper and consider this question more closely.

So, let's say we split the audio signal into a set of frames of certain length and computed related MFCC-vectors. What should we do next?


CodeBook


It is logical to assume that every frame/MFCC-vector matches with some sound. Hence, all we need is a "Sound - MFCC-vector" codebook. Once we have it we'll be able to "read" our audio data and get the words!

However, the reality is not so simple. Such codebook allows us to transform the frames/MFCC-vectors into a sequence of sounds, not letters or words. Also we shouldn't forget the following:

First of all, different people pronounce the same sound in different ways (voice quality, etc). I.e. the "Sound - MFCC-vector" relationship must be implemented as "one to many". Hopefully, that isn't a big problem and it just increases complexity of our codebook structure for a little bit.

Secondly, a sound can lasts longer than one frame. E.g. "a" sound pronounced out of a word (for recording, without forced stretching) lasts for 300ms. For consonants, this value is in 50% lower, but this doesn't change the idea.
As a result, taking into account that the frame length is 50ms, the incoming set of sounds for "odin" word could look like "ooooddiiiinn".

Thirdly, some letters might be pronounced in different ways. E.g. the english and american variant of "can't". In our example the letter "o" sometimes can be pronounced as "a" (the old russian tradition). Given that fact, we must be ready for the "aaaaddiiiinn" sequence as an incoming data.

In other words, using codebook (see "CODEBOOK" section) isn't able to completely solve the problem of speech recognition. However, it allows to reduce it to another problem, which could be perfectly solved with such a great instrument as Hidden Markov Models are.


Hidden Markov Models


There are thousands of articles [1,2] and gygabytes of videos (my favorite one) about HMM. That's why I'm not going to describe the idea in detail, but I'll say a few words about using HMM for solving the speech recognition problem. So, HMM is a set of:
  • inner/hidden states (letters in the word);
  • observed values (sounds that we get as an input data, MFCC-vectors matched by the Codebook);
  • initial distribution vector (defines from which letter the word begins);
  • transition matrix (defines order of letters in the word, how much the sound "stretches");
  • emission matrix (defines the probability to hear Y sound for X letter);
As an example let's consider word "odin" (russian "one"):
STATES  4 o d i n
OBSERVATIONS 5 o a d i n
INITIAL
1.0 0.0 0.0 0.0
TRANSITION
0.75 0.25 0.0 0.0
0.0 0.5 0.5 0.0
0.0 0.0 0.75 0.25
0.0 0.0 0.0 1.0
EMISSION
0.5 0.5 0.0 0.0 0.0
0.0 0.0 1.0 0.0 0.0
0.0 0.0 0.0 1.0 0.0
0.0 0.0 0.0 0.0 1.0

The model above is quite simple:
  • the initial distribution vector says us that word "odin" must be started from letter "o" (quite logical, eh?);
  • the transition matrix insists that vowels "o" and "i" drag on much longer, than consonants "d" and "n" (i.e. the probability of "o" to "o" transition is 0.75, in the same time "o" to "d" transition has only 0.25 probability points;
  • the emissions matrix shows us that the only variant reading is possible for "o" letter (it can be heard as "o" or "a");

Great! But what should we do with this model? Let's say we have the model above (M) and some sequence of observed sounds 0="aaaddiiinn". There are 3 standard HMM usage pattern which make it extremely useful. At the moment we need only 2 of them, so let's consider them below.


Recognition


Suppose we have a set of prepared models (a set of known words). In this case the problem of O-sequence recognition narrows to searching of the most probable model for the O-sequence. In turn, the probability of an observation sequence O given a HMM model M - P(O|M) is the probability of a sequence of observations for each possible sequence of states of the model. Yes, it sounds a little bit complicated, but, in fact, it is a simple walk through the O-sequence values (mapped by the emission matrix) and the transition matrix. The idea of this is fully described in Forward Algorithm (specification / implementation).


Learning


Now imagine that we have a set of pronunciation samples of "odin" word generated by the same person. Using this information we can fit M model parameters to increase the P(O|M) probability for each sample. In other words - we can "teach" the model and improve recognition quality. This "learning" can be done by Baum-Welch algorithm (specification / implementation).


Read more...

Friday, November 21, 2014

Использование скрытых марковских моделей для распознования речи

В предыдущей статье мы рассмотрели простейший вариант приложения для распознавания речи. Теперь я предлагаю окунуться в эту тему немного глубже и рассмотреть несколько более серьёзный подход к этом вопросу.
И так, допустим, что мы разбили входной сигнал на фреймы определённой (50мс) длины и вычислили для каждого фрейма соответствующий MFCC-вектор. Что же делать дальше?


Букварь


Логично предположить, что каждый фрейм/MFCC-вектор соответствует некоторому звуку. Следовательно, всё, что нам нужно - это составить букварь типа "Звук - MFCC-вектор". После чего мы без труда сможем "прочитать" наши данные и получить слова!

Однако, на деле всё не так просто. Подобный букварь позволит на преобразовать последовательность фреймов в последовательность звуков, а не слов. При этом не стоит забывать, что:

Во-первых, разные люди могут произносить звуки по разному: тембр голоса и т.д. То есть, отношение "Звук - MFCC-вектор" должно рассматриваться как "1 ко многим". К счастью, это небольшая проблема, которая лишь немного усложнит структуру нашего букваря.

Во-вторых, звук может длиться более одного фрейма. Так, например, краткий звук "а", произнесённый вне слова (под запись, без принудительного затягивания) длится где-то 300мс. Для согласных это значение будет где-то на 50% ниже (важно помнить, что Т - это "т", а не "тэ"), но сути дела это не меняет.
В итоге, c учётом длины фрейма в 50мс, для слова "один", мы можем получить последовательность звуков "ооооддииииннн".

В-третьих, некоторые буквы могут произноситься по разному. Так, например, буква "о", часто произносится как "а" (так называемое "аканье"). Учитывая этот факт мы вполне можем ждать на входе последовательность звуков типа "аааддииинн".

Другими словами, использование букваря (смотреть секцию "CODEBOOK") не может в полной мере решить задачу распознавания речи, однако позволяет полностью свести её в область, прекрасно подходящую для другого инструмента - Скрытых Марковских Моделей.


СММ


Про СММ написаны горы статей [1,2] и сняты гигабайты видео (моё любимое), поэтому я не буду останавливаться на них подробно, а лишь опишу способ их применения для задачи распознавания речи.

И так, скрытая марковская модель - это набор из:
  • множества внутренних/скрытых состояний (буквы в слове);
  • множества наблюдаемых значений (звуки, что мы получаем на вход из MFCC-векторов);
  • вектора начального распределения (с какой из множества заданных букв вероятнее всего начнётся слово);
  • матрицы переходов между внутренними состояниями (какая буква идёт за какой, как долго буква может "тянуться");
  • матрица эмиссий (вероятность услышать звук Y для буквы X);
В качестве примера рассмотрим слово "один":
STATES  4 o d i n
OBSERVATIONS 5 o a d i n
INITIAL
1.0 0.0 0.0 0.0
TRANSITION
0.75 0.25 0.0 0.0
0.0 0.5 0.5 0.0
0.0 0.0 0.75 0.25
0.0 0.0 0.0 1.0
EMISSION
0.5 0.5 0.0 0.0 0.0
0.0 0.0 1.0 0.0 0.0
0.0 0.0 0.0 1.0 0.0
0.0 0.0 0.0 0.0 1.0

Данная модель довольно проста:
  • вектор начального распределения говорит нам, что слово "один" всегда должно начинаться с буквы "о" (что крайне логично);
  • матрица переходов показывает, что гласные "о" и "и" тянуться намного дольше, согласных "д" и "н" (т.е. вероятность перейти из состояния "о" в него же оставляет 0.75, в то время как переход в состояние "д" имеет вероятность лишь 0.25);
  • матрица эмиссий говорит нам, что единственное разночтение возможно только для буквы "о" (может слышаться либо "о", либо "а");

Прекрасно, скажете вы, но что нам делать с этой моделью? Допустим у нас есть вышеописанная модель M и некая последовательность наблюдаемых звуков О="аааддииинн". Существует 3 стандартные задачи для СММ, 2 из которых мы сейчас расммотрим.


Распознавание


Допустим, у нас есть некоторе множество заранее подготовленных моделей (набор известных нам слов). В данном случае задача распознавания последовательности O сводится к поиску наиболее правдоподобной модели. В свою очередь правдоподобность модели для некоторой последовательности вычисляется как P(O|M) - вероятность появления последовательности наблюдений для каждой возможной последовательности состояний модели. Звучит несколько сложновато, но по факту реализуется простым "натягиванием" наблюдений на модель и оценкой того, хорошо ли "сидит". Делается это с помощью алгоритма Прямого Хода (описание / реализация).


Обучение


Теперь представим, что у нас есть несколько образцов произношения слова "один", сгенерированных конкретным человеком. Подобрав значения модели M так, что бы для всех этих последовательностей звуков модель выдавала максимальную вероятность - обучив модель максимизировав P(O|M), мы сможем значительно улучшить качество распознования. Сделать это можно с помощью алгоритма Баума-Велша (описание / реализация).


Read more...

Wednesday, June 25, 2014

Speech recognition for dummies

This article is devoted to understanding such interesting area of software development as Speech Recognition. Unfortunately, I'm not an expert in this stuff, so my clause will have a lot of inaccuracies, mistakes and disappointments. Nonetheless, the chief objective of this paper (as you can see from its topic) is not a professional analysis of the problem, but analysis of the basic concepts, issues and theirs solutions.


Prologue

Let's start from the point that our speech is a sequence of sounds. A sound, in turn, is a superposition (overlay) of acoustic oscillations (waves) of different frequencies. Finally, a wave, as we all know from physics, has two major characteristics - amplitude and frequency.

Storing a sound signal on a digital device (i.e. your computer) requires splitting the signal into a set of intervals and calculating average amplitude value for each interval.

That is how mechanical vibrations can be transformed to a set of numbers, which can be processed on our phones/computers.

Hence, the speech recognition task is reduced to "matching" a set of numerical values and words from a dictionary (i.e. Russian dictionary :).

Let's see in more details on the way of implementation of such "matching".



Input data

Suppose we have a file/stream with an audio data. First of all we need to understand its internal structure and how it can be read. Let's consider the most simple case - a WAV file.

This format implies the presence of two blocks in the file. The first one is a header with the audio stream information: bitrate, frequency, amount of channels, length of the file, etc. The second one contains the "raw" data - the digital signal, the set of signal's amplitudes values.

The reading logic in this case is quite simple. Read the header, check some restrictions (e.g. no compression), store the data into a special array.

Example.



Recognition

Theoretically, we already know all that we need to compare (element-wise) our sample (which we've read in the previous section) with another one, which text equivalent is known. In other words, we can try to recognize... But we shouldn't do this :)

Our approach must be stable (at least, a little bit) to any changes in the voice (of a man who uttering the word), volume and speed of pronunciation. Of course, we can't guarantee that with element-wise comparison of two audio signals.

That's why we choose the another way.



Frames

First of all let's split our data into a set of little time intervals – frames. Notice that frames shouldn't go one by one, instead they should “overlap” each other. Thus the end of the n-frame must intersect with the beginning of n+1-frame.

Frames are more appropriate units for analysis rather than values of singnal's amplitudes just because we are dealing with waves, and waves are supposed to be analyzed on intervals, not in particular points :) Meanwhile overlapping of frames allows us to smooth analysis results by using a frame as a some kind of “window”, which moves along the original signal (its amplitudes values).

Experiments show us that the optimal frame length should correspond to the interval of 10ms, "overlap" - to 50%. Given that the average word length is 500ms (at least in my experiments) - this length will give us approximately 500 / (10 * 0.5) = 100 frames per word.



Splitting words

The first problem, which must be solved in the process of speech recognition is splitting the speech into different words. Let's simplify our task and put that the speech contains some pauses (intervals of silence), that can be used as the words splitter.

In this case we need to find the threshold. All values above it will be words, below – silence. Here we can use several ways:

  • use a constant (the original signal must be always generated under the same circumstances, in the same way);
  • cluster signal's values into two sets: one with words values and another one with silence (works only if silence takes a big part in the original signal);
  • analyze the entropy;

As you can guess, we are going to discuss the last point :) Let's start from the entropy definition - “a measure of disorder” (c). In our case, the entropy shows how much our signal "varies" within a given frame.

In order to calculate the entropy of a given frame, perform the following steps:

  • assume that our signal is normalized and all its values ​​range in [-1;1];
  • build the histogram of the signal's values:
    • split the range into i=0..N intervals (N is a preselected value, recommended from 50 to 100);
    • compute P[i], amount of signal values, matched to the i-th interval;
    • normalize P[i] (divide on the frame's size);
  • calculate the entropy as ;

Example.

So, we got the value of the entropy. But this is just an another frame's characteristic. We still need to compare it with something if we want to separate the sound and the silence.
Fortunately, entropy (in contrast to Root Mean Square) is a relatively independent measure. This fact allowed us to put its threashold value as a constant (e.g. 0.1).

Nevertheless, the problems don't end here :( The entropy may sag in the middle of a word (in vowels) or may suddenly jump (in case of a little noise). In order to fix the first problem we need to use the “minimal distance between two words” and “glue” nearby frame sets, which were divided because of the sag. The second problem can be resolved by using “minimal word length” and removing all the candidates which are not as long as we need (and aren't used in the first point).

If the speech doesn't have any pauses between words (or the speaker speaks too fast), we can try to create/extract sets of various subsequences from the original frame set. Then all these subsequences can be used as a source for the recognition process. But that's absolutely another story :)



MFCC

So, we have the set of frames, which matches with a particular word. We can take the path of least resistance and use average square (root mean square) of all the frame's elements as theirs numerical characteristic. However, such metric gives us a quite little information for further analysis.

That's why it's time to introduce Mel-Frequency Cepstral Coefficients. According to Wikipedia (which, as we know, never lies) MFCC are a kind of signal's energy representation. Pros of their usage are as follows:

  • using of the spectrum of the signal (i.e. expansion by the orthogonal basis of [co]sine functions), which takes into account the “wave” nature of the signal;
  • the spectrum is projected onto a special mel-scale, which allows us to extract the most valuable for human perception frequencies;
  • amount of calculated coefficients may be limited by any value (e.g. 12), that allows us to “squeeze” the frame and, as a consequence, amount of information to further processing;

It's time to consider the process of MFCC computing for a particular frame.

Let's think about our frame as a vector , where N is its size.



Fourier transformation

First of all let's calculate the spectrum of the signal by computing the discreet Fourier transformation (preferably its “fast” FFT implementation).

It's also recommended to apply so-called “window” Hamming function to “smooth” values on the frame's bounds.

As a result we will have the following vector:

It is important to understand that after this conversion, we have X-axis for the frequency (hz) signal, and the Y-axis for the magnitude (as a way to escape from the complex values​​):



Mel-coefficients calculation

Let's start from the “mel” definition. According to Wikipedia (again), mel is a "psychophysical unit of sound pitch, based on the subjective perception of an average human. It primarily depends on the frequency of the sound (as well as on volume and tone). In other words, this value shows us how much the sound of a certain frequency is "valuable" for us.

Frequency can be converted into mel-values by the following formula (remind it as "formula-1"):

Backward transformation looks in the following way (remind it as "formula-2"):

Graph of mel/frequency scale:

But let's get back to our task. Assume that we have a frame which size is N = 256 elements. We know (from the audio format data) that the frequency of the frame is 16000hz. Let's put that the human speech belongs to the range [300; 8000]hz. Amount of mel-coefficients that we want to find is M = 10.

In order to expand the spectrum (which we've calculated above) by the mel-scale, we need to create the “comb” of mel-filters. In fact, each mel-filter is a triangular window function, which summarizes amount of energy at the certain frequency range (thus calculates the mel-coefficient). Hence, since we know the amount of mel-coefficients and the frequency range, we can construct the following set of filters:

Notice, the greater the index number of the mel-coefficient is, the wider the base of the filter is. This is because the process of splitting frequency bands into the filters intervals happens on the mel-scale (not on the frequency-scale).

But I digress. For our case the frequency range is [300, 8000]. According to the “formula-1” on the mel-scale this interval transforms into [401.25; 2834.99].

Next, we need 12 reference points to build 10 triangular filters:

m[i] = [401.25, 622.50, 843.75, 1065.00, 1286.25, 1507.50, 1728.74, 1949.99, 2171.24, 2392.49, 2613.74, 2834.99]

Notice, on the mel-scale the points are located uniformly. Let's transform the values back to frequency using the “formula-2”:

h[i] = [300, 517.33, 781.90, 1103.97, 1496.04, 1973.32, 2554.33, 3261.62, 4122.63, 5170.76, 6446.70, 8000]

As you can see, our values start to stretch, thus aligning the growth dynamics of the "significance" of the sound on the low and high frequencies.

Now let's impose the scale (which we got above) on the spectrum of our frame. As we remember, X-axis is for frequency. Length of the spectrum is 256 elements, which refers to 16000Hz diapason. Solving the simple proportion we obtain the following formula:

f(i) = floor( (frameSize+1) * h(i) / sampleRate)

in our case this is equivalent to

f(i) = 4, 8, 12, 17, 23, 31, 40, 52, 66, 82, 103, 128

That is all! Since we know the reference points on the X-axis of our spectrum we can easily construct the necessary filters by the following formula:



Applying the mel-filters

Applying the filter to the spectrum mean pairwise multiplying the filter values ​​with the spectrum values. The sum of the elements in the result vector is a mel-coefficient. Since we have M-filters, the number of mel-coefficients will be the same.

However, we need to apply mel-filters not to the spectrum, but to its energy. Then we must logarithm the result. It's believed that the trick above decreases the sensitivity to noise ratios.



Cosine transformation

Discrete Cosine Transform (DCT) is used in order to get the "cepstral" coefficients. It allows to “squeeze” the results, increasing the importance of the first coefficients and reducing the importance of the last ones.

In current case we use DCT-II without additional multiplication on the scale factor .

Now we have the set of M mfcc-coefficients for each frame. These coefficients can (and will) be used for the further analysis!

Examples of the methods above can be found here.



Recognition algorithm

Here, dear reader, you'll face the major disappointment. In the internets I've seen a lot of highly (and not very) intelligent debates about what is the best way for speech recognizing. Somebody stands up for Hidden Markov Models, someone - for neural networks... someone's thoughts can't be understood at all :)

In any case a lot of people choose HMM, and this is the way which I'm going to implement… in the future :)

At the moment I suggest to stop on much less effective, but much more simple way.

So, let's remember that our task is about recognition words from a dictionary. For simplicity's sake, we will recognize first ten numbers from Russian dictionary: “odin”, “dva”, “tri”, “chetyre”, “pyat”, “shest”, “sem”, “vosem”, “devyat”, “desyat” :)

Now, take your iPhone/Android and walk through your colleagues asking them to dictate these words for recording. Save the audio files on your computer and convert them into Wav format. Next, calculate mfcc-coefficients for each audio file (a set/vector of concatenated mfcc-coefficients of all frames from the file) and store matches between words and related sets of mfcc-coefficients into a local storage (file or database);

Let's call each [word; sets of mfcc-coefficients] match “Model”, and the process of creating these matches - “Machine Learning”! In fact, the simple adding of new samples in the local storage has a very tenuous connection with machine learning ... But it's a painfully trendy term :)

Now our task is reduced to the selection of the most "close" model for the set of mfcc-coefficients (the word which we want to recognize). At first glance, the problem can be solved quite simply:

  • calculate distances between each model and the word we want to recognize (distance between the model and the word is average euclid distance between word's mfcc-coefficients vector and model's mfcc-coefficients vectors);
  • choose as a correct the model with the shortest distance;

However, two different speakers can pronounce the same word with different speed. That means size of mfcc-vector can be different for two sample of the same word.

Fortunately, the problem of comparing sequences of different length has been solved in the form of Dynamic Time Warping algorithm. This dynamic programming algorithm is perfectly described on Wiki.

The only thing that we need to do with that algorithm is changing of the distance calculation logic. We must remember that the mfcc-vector of the model is a sequence of mfcc-subvectors (M-size) obtained from the frames. So, DTW algorithm must calculate distances between these subvectors, not their elements. Hence, elements of distance matrix are Euclid distances between frame based mfcc-subvectors.

Example.



Experiments

I haven't been able to verify the approach on a big amount of “training” samples. As for test results based on 3 samples education... They are quite bad – only 65% correct recognitions.

However, my task was implementation of the simplest speech recognition application. So called “proof of concept” :)



Implementation

The attentive reader has noticed that the article contains a lot of references on the GitHub project. It's worth to say that it's my first C++ project since the university times. Also it's my first attempt to calculate something more complicated than the arithmetic mean since the same times... In other words “it comes with absolutely no warranty” (с) :)




Read more...

Friday, June 13, 2014

Распознование речи для чайников

В этой статье я хочу рассмотреть основы такой интереснейшей области разработки ПО как Распознавание Речи. Экспертом в данной теме я, естественно, не являюсь, поэтому мой рассказ будет изобиловать неточностями, ошибками и разочарованиями. Тем не менее, главной целью моего "труда", как можно понять из названия, является не профессиональный разбор проблемы, а описание базовых понятий, проблем и их решений. В общем, прошу всех заинтересовавшихся пожаловать под кат!



Пролог

Начнём с того, что наша речь - это последовательность звуков. Звук в свою очередь - это суперпозиция (наложение) звуковых колебаний (волн) различных частот. Волна же, как нам известно из физики, характеризуются двумя атрибутами - амплитудой и частотой.

Для того, что бы сохранить звуковой сигнал на цифровом носителе, его необходимо разбить на множество промежутков и взять некоторое "усредненное" значение на каждом из них.

Таким вот образом механические колебания превращаются в набор чисел, пригодный для обработки на современных ЭВМ.

Отсюда следует, что задача распознавания речи сводится к "сопоставлению" множества численных значений (цифрового сигнала) и слов из некоторого словаря (русского языка, например).

Давайте разберемся, как, собственно, это самое "сопоставление" может быть реализовано.



Входные данные

Допустим у нас есть некоторый файл/поток с аудиоданными. Прежде всего нам нужно понять, как он устроен и как его прочесть. Давайте рассмотрим самый простой вариант - WAV файл.

Формат подразумевает наличие в файле двух блоков. Первый блок - это заголовка с информацией об аудиопотоке: битрейте, частоте, количестве каналов, длине файла и т.д. Второй блок состоит из "сырых" данных - того самого цифрового сигнала, набора значений амплитуд.

Логика чтения данных в этом случае довольно проста. Считываем заголовок, проверяем некоторые ограничения (отсутствие сжатия, например), сохраняем данные в специально выделенный массив.

Пример.



Распознавание

Чисто теоретически, теперь мы можем сравнить (поэлементно) имеющийся у нас образец с каким-нибудь другим, текст которого нам уже известен. То есть попробовать "распознать" речь... Но лучше этого не делать :)

Наш подход должен быть устойчив (ну хотя бы чуть-чуть) к изменению тембра голоса (человека, произносящего слово), громкости и скорости произношения. Поэлементным сравнением двух аудиосигналов этого, естественно, добиться нельзя.

Поэтому мы пойдем несколько иным путём.



Фреймы

Первым делом разобьём наши данные по небольшим временным промежуткам - фреймам. Причём фреймы должны идти не строго друг за другом, а “внахлёст”. Т.е. конец одного фрейма должен пересекаться с началом другого.

Фреймы являются более подходящей единицей анализа данных, чем конкретные значения сигнала, так как анализировать волны намного удобней на некотором промежутке, чем в конкретных точках. Расположение же фреймов “внахлёст” позволяет сгладить результаты анализа фреймов, превращая идею фреймов в некоторое “окно”, движущееся вдоль исходной функции (значений сигнала).

Опытным путём установлено, что оптимальная длина фрейма должна соответствовать промежутку в 10мс, "нахлёст" - 50%. С учётом того, что средняя длина слова (по крайней мере в моих экспериментах) составляет 500мс - такой шаг даст нам примерно 500 / (10 * 0.5) = 100 фреймов на слово.



Разбиение слов

Первой задачей, которую приходится решать при распознавании речи, является разбиение этой самой речи на отдельные слова. Для простоты предположим, что в нашем случае речь содержит в себе некоторые паузы (промежутки тишины), которые можно считать “разделителями” слов.

В таком случае нам нужно найти некоторое значение, порог - значения выше которого являются словом, ниже - тишиной. Вариантов тут может быть несколько:

  • задать константой (сработает, если исходный сигнал всегда генерируется при одних и тех же условиях, одним и тем же способом);
  • кластеризовать значения сигнала, явно выделив множество значений соответствующих тишине (сработает только если тишина занимает значительную часть исходного сигнала);
  • проанализировать энтропию;

Как вы уже догадались, речь сейчас пойдёт о последнем пункте :) Начнём с того, что энтропия - это мера беспорядка, “мера неопределённости какого-либо опыта” (с). В нашем случае энтропия означает то, как сильно “колеблется” наш сигнал в рамках заданного фрейма.

Для того, что бы подсчитать энтропию конкретного фрейма выполним следующие действия:

  • предположим, что наш сигнал пронормирован и все его значения лежат в диапазоне [-1;1];
  • построим гистограмму (плотность распределения) значений сигнала фрейма:
    • разобьём этот диапазон на i=0..N промежутков (N-выбирается заранее, рекомендуется от 50 до 100);
    • посчитаем P[i]-ые, количество значений сигнала, входящих в i-ый промежуток;
    • пронормируем P[i] (разделим на размер фрейма);
  • рассчитаем энтропию, как ;

Пример.

И так, мы получили значение энтропии. Но это всего лишь ещё одна характеристика фрейма, и для того, что бы отделить звук от тишины, нам по прежнему нужно её с чем-то сравнивать. В некоторых статьях рекомендуют брать порог энтропии равным среднему между её максимальным и минимальным значениями (среди всех фреймов). Однако, в моём случае такой подход не дал сколь либо хороших результатов.
К счастью, энтропия (в отличие от того же среднего квадрата значений) - величина относительно самостоятельная. Что позволило мне подобрать значение её порога в виде константы (0.1).

Тем не менее проблемы на этом не заканчиваются :( Энтропия может проседать по середине слова (на гласных), а может внезапно вскакивать из-за небольшого шума. Для того, что бы бороться с первой проблемой, приходится вводить понятие “минимально расстояния между словами” и “склеивать” близ лежачие наборы фреймов, разделённые из-за проседания. Вторая проблема решается использованием “минимальной длины слова” и отсечением всех кандидатов, не прошедших отбор (и не использованных в первом пункте).

Если же речь в принципе не является “членораздельной”, можно попробовать разбить исходный набор фреймов на определённым образом подготовленные подпоследовательности, каждая из которых будет подвергнута процедуре распознавания. Но это уже совсем другая история :)



MFCC

И так, мы у нас есть набор фреймов, соответствующих определённому слову. Мы можем пойти по пути наименьшего сопротивления и в качестве численной характеристики фрейма использовать средний квадрат всех его значений (Root Mean Square). Однако, такая метрика несёт в себе крайне мало пригодной для дальнейшего анализа информации.

Вот тут в игру и вступают Мел-частотные кепстральные коэффициенты (Mel-frequency cepstral coefficients). Согласно Википедии (которая, как известно, не врёт) MFCC - это своеобразное представление энергии спектра сигнала. Плюсы его использования заключаются в следующем:

  • Используется спектр сигнала (то есть разложение по базису ортогональных [ко]синусоидальных функций), что позволяет учитывать волновую “природу” сигнала при дальнейшем анализе;
  • Спектр проецируется на специальную mel-шкалу, позволяя выделить наиболее значимые для восприятия человеком частоты;
  • Количество вычисляемых коэффициентов может быть ограничено любым значением (например, 12), что позволяет “сжать” фрейм и, как следствие, количество обрабатываемой информации;

Давайте рассмотрим процесс вычисления MFCC коэффициентов для некоторого фрейма.

Представим наш фрейм в виде вектора , где N - размер фрейма.



Разложение в ряд Фурье

Первым делом рассчитываем спектр сигнала с помощью дискретного преобразования Фурье (желательно его “быстрой” FFT реализацией).

Так же к полученным значениям рекомендуется применить оконную функцию Хэмминга, что бы “сгладить” значения на границах фреймов.

То есть результатом будет вектор следующего вида:

Важно понимать, что после этого преобразования по оси Х мы имеем частоту (hz) сигнала, а по оси Y - магнитуду (как способ уйти от комплексных значений):



Расчёт mel-фильтров

Начнём с того, что такое mel. Опять же согласно Википедии, mel - это “психофизическая единица высоты звука”, основанная на субъективном восприятии среднестатистическими людьми. Зависит в первую очередь от частоты звука (а так же от громкости и тембра). Другими словами, эта величина, показывающая, на сколько звук определённой частоты “значим” для нас.

Преобразовать частоту в мел можно по следующей формуле (запомним её как "формула-1"):

Обратное преобразование выглядит так (запомним её как "формула-2"):

График зависимости mel / частота:

Но вернёмся к нашей задаче. Допустим у нас есть фрейм размером 256 элементов. Мы знаем (из данных об аудиоформате), что частота звука в данной фрейме 16000hz. Предположим, что человеческая речь лежит в диапазоне от [300; 8000]hz. Количество искомых мел-коэффициентов положим M = 10 (рекомендуемое значение).

Для того, что бы разложить полученный выше спектр по mel-шкале, нам потребуется создать “гребёнку” фильтров. По сути, каждый mel-фильтр это треугольная оконная функция, которая позволяет просуммировать количество энергии на определённом диапазоне частот и тем самым получить mel-коэффициент. Зная количество мел-коэффициентов и анализируемый диапазон частот мы можем построить набор таких вот фильтров:

Обратите внимание, что чем больше порядковый номер мел-коэффициента, тем шире основание фильтра. Это связано с тем, что разбиение интересующего нас диапазона частот на обрабатываемые фильтрами диапазоны происходит на шкале мелов.

Но мы опять отвлеклись. И так для нашего случая диапазон интересующих нас частот равен [300, 8000]. Согласно формуле-1 в на мел-шкале этот диапазон превращается в [401.25; 2834.99].

Далее, для того, что бы построить 10 треугольных фильтров нам потребуется 12 опорных точек:

m[i] = [401.25, 622.50, 843.75, 1065.00, 1286.25, 1507.50, 1728.74, 1949.99, 2171.24, 2392.49, 2613.74, 2834.99]

Обратите внимание, что на мел-шкале точки расположены равномерно. Переведём шкалу обратно в герцы с помощью формулы-2:

h[i] = [300, 517.33, 781.90, 1103.97, 1496.04, 1973.32, 2554.33, 3261.62, 4122.63, 5170.76, 6446.70, 8000]

Как видите теперь шкала стала постепенно растягиваться, выравнивая тем самым динамику роста “значимости” на низких и высоких частотах.

Теперь нам нужно наложить полученную шкалу на спектр нашего фрейма. Как мы помним, по оси Х у нас находится частота. Длина спектра 256 - элементов, при этом в него умещается 16000hz. Решив нехитрую пропорцию можно получить следующую формулу:

f(i) = floor( (frameSize+1) * h(i) / sampleRate)

что в нашем случае эквивалентно

f(i) = 4, 8, 12, 17, 23, 31, 40, 52, 66, 82, 103, 128

Вот и всё! Зная опорные точки на оси Х нашего спектра, легко построить необходимые нам фильтры по следующей формуле:



Применение фильтров, логарифмирование энергии спектра

Применение фильтра заключается в попарном перемножении его значений со значениями спектра. Результатом этой операции является mel-коэффициент. Поскольку фильтров у нас M, коэффициентов будет столько же.

Однако, нам нужно применить mel-фильтры не к значениям спектра, а к его энергии. После чего прологарифмировать полученные результаты. Считается, что таким образом понижается чувствительность коэффициентов к шумам.



Косинусное преобразование

Дискретное косинусное преобразование (DCT) используется для того, что бы получить те самые “кепстральные” коэффициенты. Смысл его в том, что бы “сжать” полученные результаты, повысив значимость первых коэффициентов и уменьшив значимость последних.

В данном случае используется DCTII без каких-либо домножений на (scale factor).

Теперь для каждого фрейма мы имеем набор из M mfcc-коэффициентов, которые могут быть использованы для дальнейшего анализа.

Примеры код для вышележащих методов можно найти тут.



Алгоритм распознавания

Вот тут, дорогой читатель, тебя и ждёт главное разочарование. В интернетах мне довелось увидеть множество высокоинтеллектуальных (и не очень) споров о том, какой же способ распознавания лучше. Кто-то ратует за Скрытые Марковские Модели, кто-то - за нейронные сети, чьи-то мысли в принципе невозможно понять :)

В любом случае немало предпочтений отдаётся именно СММ, и именно их реализацию я собираюсь добавить в свой код… в будущем :)

На данный момент, предлагаю остановится на гораздо менее эффективном, но в разы более простом способе.

И так, вспомним, что наша задача заключается в распознавании слова из некоторого словаря. Для простоты, будем распознавать называния первых десять цифр: “один“, “два“, “три“, “четыре“, “пять“, “шесть“, “семь“, “восемь“, “девять“, “десять“.

Теперь возьмем в руки айфон/андроид и пройдёмся по L коллегам с просьбой продиктовать эти слова под запись. Далее поставим в соответствие (в какой-нибудь локальной БД или простом файле) каждому слову L наборов mfcc-коэффициентов соответствующих записей.

Это соответствие мы назовём “Модель”, а сам процесс - Machine Learning! На самом деле простое добавление новых образцов в базу имеет крайне слабую связь с машинным обучением… Но уж больно термин модный :)

Теперь наша задача сводится к подбору наиболее “близкой” модели для некоторого набора mfcc-коэффициентов (распознаваемого слова). На первый взгляд задачу можно решить довольно просто:

  • для каждой модели находим среднее (евклидово) расстояние между идентифицируемым mfcc-вектором и векторами модели;
  • выбираем в качестве верной ту модель, среднее расстояние до которой будет наименьшим;

Однако, одно и тоже слово может произносится как Андреем Малаховым, так и каким-нибудь его эстонским коллегой. Другими словами размер mfcc-вектора для одного и того же слова может быть разный.

К счастью, задача сравнения последовательностей разной длины уже решена в виде Dynamic Time Warping алгоритма. Этот алгоритм динамическо программирования прекрасно расписан как в буржуйской Wiki, так и на православном Хабре.

Единственное изменение, которое в него стоит внести - это способ нахождения дистанции. Мы должны помнить, что mfcc-вектор модели - на самом деле последовательность mfcc-“подвекторов” размерности M, полученных из фреймов. Так вот, DTW алгоритм должен находить дистанцию между последовательностями эти самых “подвекторов” размерности M. То есть в качестве значений матрицы расстояний должны использовать расстояния (евклидовы) между mfcc-“подвекторами” фреймов.

Пример.



Эксперименты

У меня не было возможности проверить работу данного подхода на большой “обучающей” выборке. Результаты же тестов на выборке из 3х экземпляров для каждого слова в несинтетических условиях показали мягко говоря нелучший результат - 65% верных распознаваний.

Тем не менее моей задачей было создание максимального простого приложения для распознавания речи. Так сказать “proof of concept” :)



Реализация

Внимательный читатель заметил, что статья содержит множество ссылок на GitHub-проект. Тут стоит отметить, что это мой первый проект на С++ со времён университета. Так же это моя первая попытка рассчитать что-то сложнее среднего арифметического со времён того же университета... Другими словами it comes with absolutely no warranty (с) :)




Read more...